knitr::opts_chunk$set(echo = TRUE)

rm(list = ls())

#cmdstanr::set_cmdstan_path(path = "C:/Users/kueng/.cmdstan/cmdstan-2.35.0")
#cmdstanr::set_cmdstan_path(path = "C:/Users/pascku/.cmdstan/cmdstan-2.36.0")

library(tidyverse)
library(R.utils)
library(wbCorr)
library(readxl)
library(kableExtra)
library(brms)
library(bayesplot)
library(see)
library(beepr)
library(DHARMa)
library(digest)



source(file.path('Functions', 'ReportModels.R'))
source(file.path('Functions', 'PrettyTables.R'))
source(file.path('Functions', 'ReportMeasures.R'))
source(file.path('Functions', 'PrepareData.R'))
system("shutdown /a")
## [1] 1116
# Set options for analysis
use_mi = FALSE
shutdown = FALSE
report_ordinal = FALSE
do_priorsense = FALSE
get_bayesfactor = TRUE
check_models = TRUE #

if (get_bayesfactor) {
  stats_to_report <- c('CI', 'SE', 'pd', 'ROPE', 'BF', 'Rhat', 'ESS')
} else {
  stats_to_report <- c('CI', 'SE', 'pd', 'ROPE', 'Rhat', 'ESS')
}

options(
  dplyr.print_max = 100, 
  brms.backend = 'cmdstan',
  brms.file_refit = ifelse(use_mi, 'never', 'on_change'),
  brms.file_refit = 'on_change',
  #brms.file_refit = 'always',
  error = function() {
    beepr::beep(sound = 5)
    if (shutdown) {
      system("shutdown /s /t 180")
      quit(save = "no", status = 1)
    }
  }
  , es.use_symbols = TRUE
)


####################### Model parameters #######################

iterations = 2000 # 12'000 per chain to achieve 40'000
warmup = 1000 # 2000

# NO AR!!!
#corstr = 'ar'
#corstr = 'cosy_couple'
#corstr = 'cosy_couple:user'


################################################################

suffix = paste0('_as_preregistered_new_', as.character(iterations))
df <- openxlsx::read.xlsx(file.path('long.xlsx'))
df_original <- df

df_double <- prepare_data(df, recode_pushing = TRUE, use_mi = use_mi)[[1]]

Constructing scales Re-coding pusing reshaping data (4field) centering data within and between

summary(df_double$pushing)

Min. 1st Qu. Median Mean 3rd Qu. Max. NA’s 0.0000 0.0000 0.0000 0.1649 0.0000 5.0000 275

Prepare for APIM

df_double$is_A <- ifelse(df_double$pa_obj_self_cb >= df_double$pa_obj_partner_cb, 0, 1)
df_double$is_B <- ifelse(df_double$pa_obj_self_cb < df_double$pa_obj_partner_cb, 0, 1)

Modelling

Self-Reported MVPA

Persuasion

formula <- bf(
  pa_sub ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:persuasion_self_cw + is_B:persuasion_self_cw + 
    is_A:persuasion_partner_cw + is_B:persuasion_partner_cw +
    
    is_A:persuasion_self_cb + is_B:persuasion_self_cb + 
    is_A:persuasion_partner_cb + is_B:persuasion_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    is_A:day + is_B:day + 
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:persuasion_self_cw + is_B:persuasion_self_cw + 
    is_A:persuasion_partner_cw + is_B:persuasion_partner_cw |dd| coupleID),
  
  hu = ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:persuasion_self_cw + is_B:persuasion_self_cw + 
    is_A:persuasion_partner_cw + is_B:persuasion_partner_cw +
    
    is_A:persuasion_self_cb + is_B:persuasion_self_cb + 
    is_A:persuasion_partner_cb + is_B:persuasion_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    is_A:day + is_B:day + 
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:persuasion_self_cw + is_B:persuasion_self_cw + 
    is_A:persuasion_partner_cw + is_B:persuasion_partner_cw |dd| coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
) 

prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 2)", class = "b", dpar = "hu")
  
  , brms::set_prior("normal(0, 50)", class = "b", coef = 'is_A')
  , brms::set_prior("normal(0, 50)", class = "b", coef = 'is_B')
  
  , brms::set_prior("normal(0.5, 2.5)", class = "b", coef = 'is_A', dpar = 'hu')
  , brms::set_prior("normal(0.5, 2.5)", class = "b", coef = 'is_B', dpar = 'hu')
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = hurdle_lognormal()
#)

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_sub_persuasion <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::hurdle_lognormal(), 
  #family = brms::hurdle_negbinomial(), 
  #family = brms::hurdle_poisson(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 42,
  file = file.path("models_cache_brms", paste0("pa_sub_hu_lognormal_persuasion", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
if (check_models) {
  check_brms(pa_sub_persuasion, log_pp_check = TRUE)
  DHARMa.check_brms.all(pa_sub_persuasion, integer = TRUE, outliers_type = 'bootstrap')
}
## Model has no intercept. VIFs may not be sensible.
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                        Term  VIF   VIF 95% CI Increased SE Tolerance
##                        is_A 2.24 [2.15, 2.34]         1.50      0.45
##                        is_B 1.06 [1.04, 1.10]         1.03      0.94
##     is_A:persuasion_self_cw 1.17 [1.14, 1.21]         1.08      0.86
##     is_B:persuasion_self_cw 1.15 [1.12, 1.19]         1.07      0.87
##  is_A:persuasion_partner_cw 1.05 [1.03, 1.10]         1.03      0.95
##  is_B:persuasion_partner_cw 2.78 [2.66, 2.90]         1.67      0.36
##     is_A:persuasion_self_cb 2.83 [2.70, 2.96]         1.68      0.35
##  Tolerance 95% CI
##      [0.43, 0.46]
##      [0.91, 0.96]
##      [0.83, 0.88]
##      [0.84, 0.90]
##      [0.91, 0.97]
##      [0.34, 0.38]
##      [0.34, 0.37]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        cauchy         94%
##     lognormal          6%
## 
## Predicted Distribution of Response
## 
##                Distribution Probability
##  neg. binomial (zero-infl.)         78%
##               beta-binomial         16%
##                   lognormal          3%
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 6, observations = 3742, p-value = 1
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0.0005344735 0.0030799038
## sample estimates:
## outlier frequency (expected: 0.00170229823623731 ) 
##                                        0.001603421
# rope range for continuous part of the model
rope_factor <- sd(log(pa_sub_persuasion$data$pa_sub[pa_sub_persuasion$data$pa_sub > 0]))
rope_range_continuous = c(-0.1 * rope_factor, 0.1 * rope_factor)

summary_pa_sub <- summarize_brms(
  pa_sub_persuasion, 
  stats_to_report = stats_to_report,
  rope_range = rope_range_continuous,
  hu_rope_range = c(-0.18, 0.18),
  #model_rows_fixed = model_rows_fixed,
  #model_rows_random = model_rows_random,
  #model_rownames_fixed = model_rownames_fixed,
  #model_rownames_random = model_rownames_random,
  exponentiate = T) 
## Model has no intercept. VIFs may not be sensible.
## Model has no intercept. VIFs may not be sensible.
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior
##   samples is recommended.
## Warning in summarize_brms(pa_sub_persuasion, stats_to_report =
## stats_to_report, : Coefficients were exponentiated. Double check if this
## was intended.
# Print the updated dataframe
summary_pa_sub %>%
  print_df(
    #rows_to_pack = rows_to_pack
    )
exp(Est.)_hu SE_hu 95% CI_hu pd_hu ROPE_hu inside ROPE_hu BF_hu BF_Evidence_hu Rhat_hu Bulk_ESS_hu Tail_ESS_hu exp(Est.)_nonzero SE_nonzero 95% CI_nonzero pd_nonzero ROPE_nonzero inside ROPE_nonzero BF_nonzero BF_Evidence_nonzero Rhat_nonzero Bulk_ESS_nonzero Tail_ESS_nonzero
is_A 0.74 0.14 [ 0.51, 1.08] 0.939 [0.84, 1.20] 0.270 0.270 Moderate Evidence for Null 1.002 2269 2450 46.64*** 3.64 [39.95, 53.88] 1.000 [0.92, 1.08] 0.000 >100 Overwhelming Evidence 1.000 1427 2428
is_A:day 0.97 0.18 [ 0.67, 1.41] 0.575 [0.84, 1.20] 0.650 0.098 Strong Evidence for Null 1.001 5389 2918 1.09 0.09 [ 0.92, 1.29] 0.851 [0.92, 1.08] 0.434 0.011 Very Strong Evidence for Null 1.002 4676 2918
is_A:persuasion_partner_cb 0.76 0.30 [ 0.32, 1.62] 0.765 [0.84, 1.20] 0.286 0.261 Moderate Evidence for Null 1.001 1551 2445 1.03 0.17 [ 0.75, 1.42] 0.568 [0.92, 1.08] 0.365 0.017 Very Strong Evidence for Null 1.002 1518 2131
is_A:persuasion_partner_cw 1.58*** 0.16 [ 1.29, 1.97] 1.000 [0.84, 1.20] 0.004 >100 Overwhelming Evidence 1.000 2783 3038 1.02 0.03 [ 0.96, 1.08] 0.716 [0.92, 1.08] 0.981 0.008 Very Strong Evidence for Null 1.001 2593 2466
is_A:persuasion_self_cb 0.81 0.38 [ 0.33, 1.99] 0.677 [0.84, 1.20] 0.274 0.268 Moderate Evidence for Null 1.000 1243 2299 1.21 0.21 [ 0.84, 1.73] 0.866 [0.92, 1.08] 0.190 0.032 Strong Evidence for Null 1.002 1506 2250
is_A:persuasion_self_cw 1.67*** 0.16 [ 1.41, 2.09] 1.000 [0.84, 1.20] 0.000 >100 Overwhelming Evidence 1.002 2279 1744 1.01 0.04 [ 0.94, 1.09] 0.650 [0.92, 1.08] 0.935 0.011 Very Strong Evidence for Null 1.005 1848 2443
is_B 0.92 0.18 [ 0.63, 1.37] 0.657 [0.84, 1.20] 0.589 0.083 Strong Evidence for Null 1.002 2113 2747 49.88*** 3.87 [42.42, 57.98] 1.000 [0.92, 1.08] 0.000 >100 Overwhelming Evidence 1.004 1498 2350
is_B:day 0.87 0.16 [ 0.61, 1.27] 0.776 [0.84, 1.20] 0.542 0.122 Moderate Evidence for Null 1.002 5012 2585 0.88 0.07 [ 0.75, 1.03] 0.946 [0.92, 1.08] 0.246 0.027 Very Strong Evidence for Null 1.000 5060 2971
is_B:persuasion_partner_cb 1.01 0.50 [ 0.38, 2.59] 0.507 [0.84, 1.20] 0.276 0.246 Moderate Evidence for Null 1.001 1212 1923 1.03 0.20 [ 0.71, 1.48] 0.562 [0.92, 1.08] 0.319 0.016 Very Strong Evidence for Null 1.003 1293 1942
is_B:persuasion_partner_cw 1.48*** 0.14 [ 1.24, 1.85] 1.000 [0.84, 1.20] 0.011 37.706 Very Strong Evidence 1.001 3118 1930 1.03 0.03 [ 0.97, 1.10] 0.850 [0.92, 1.08] 0.929 0.015 Very Strong Evidence for Null 1.001 2987 2511
is_B:persuasion_self_cb 0.70 0.30 [ 0.30, 1.58] 0.797 [0.84, 1.20] 0.234 0.282 Moderate Evidence for Null 1.001 1614 2463 1.13 0.18 [ 0.81, 1.56] 0.769 [0.92, 1.08] 0.271 0.023 Very Strong Evidence for Null 1.002 1431 2204
is_B:persuasion_self_cw 1.68*** 0.17 [ 1.39, 2.12] 1.000 [0.84, 1.20] 0.001 >100 Overwhelming Evidence 1.000 2613 2076 1.03 0.03 [ 0.98, 1.10] 0.876 [0.92, 1.08] 0.932 0.015 Very Strong Evidence for Null 1.001 2680 3144
sd(is_A) 0.94 0.13 [0.73, 1.23] NA NA NA NA NA 1.002 1960 2814 0.33 0.05 [0.25, 0.45] NA NA NA NA NA 1.001 1331 2454
sd(is_A:persuasion_partner_cw) 0.42 0.11 [0.21, 0.70] NA NA NA NA NA 1.000 1999 2257 0.06 0.03 [0.01, 0.13] NA NA NA NA NA 1.006 1037 757
sd(is_A:persuasion_self_cw) 0.29 0.13 [0.05, 0.59] NA NA NA NA NA 1.003 1124 1201 0.17 0.03 [0.11, 0.23] NA NA NA NA NA 1.004 2044 2529
sd(is_B) 1.01 0.13 [0.79, 1.32] NA NA NA NA NA 1.001 2288 2551 0.39 0.05 [0.30, 0.50] NA NA NA NA NA 1.001 1665 2636
sd(is_B:persuasion_partner_cw) 0.33 0.14 [0.05, 0.63] NA NA NA NA NA 1.003 1118 958 0.10 0.03 [0.05, 0.17] NA NA NA NA NA 1.003 2312 2581
sd(is_B:persuasion_self_cw) 0.38 0.11 [0.19, 0.62] NA NA NA NA NA 1.000 2229 2918 0.09 0.03 [0.05, 0.16] NA NA NA NA NA 1.001 2289 2531
sigma NA NA NA NA NA NA NA NA NA NA NA 0.67 0.01 [0.65, 0.70] NA NA NA NA NA 1.001 5885 2983

Hurdle part of the model on the left, non-zero part towards the right side of the table

Pressure

formula <- bf(
  pa_sub ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:pressure_self_cw + is_B:pressure_self_cw + 
    is_A:pressure_partner_cw + is_B:pressure_partner_cw +
    
    is_A:pressure_self_cb + is_B:pressure_self_cb + 
    is_A:pressure_partner_cb + is_B:pressure_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    is_A:day + is_B:day + 
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:pressure_self_cw + is_B:pressure_self_cw + 
    is_A:pressure_partner_cw + is_B:pressure_partner_cw |dd| coupleID),
  
  hu = ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:pressure_self_cw + is_B:pressure_self_cw + 
    is_A:pressure_partner_cw + is_B:pressure_partner_cw +
    
    is_A:pressure_self_cb + is_B:pressure_self_cb + 
    is_A:pressure_partner_cb + is_B:pressure_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    is_A:day + is_B:day + 
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:pressure_self_cw + is_B:pressure_self_cw + 
    is_A:pressure_partner_cw + is_B:pressure_partner_cw |dd| coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
) 

prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 2)", class = "b", dpar = "hu")
  
  , brms::set_prior("normal(0, 50)", class = "b", coef = 'is_A')
  , brms::set_prior("normal(0, 50)", class = "b", coef = 'is_B')
  
  , brms::set_prior("normal(0.5, 2.5)", class = "b", coef = 'is_A', dpar = 'hu')
  , brms::set_prior("normal(0.5, 2.5)", class = "b", coef = 'is_B', dpar = 'hu')
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = hurdle_lognormal()
#)

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_sub_pressure <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::hurdle_lognormal(), 
  #family = brms::hurdle_negbinomial(), 
  #family = brms::hurdle_poisson(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 42,
  file = file.path("models_cache_brms", paste0("pa_sub_hu_lognormal_pressure", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
if (check_models) {
  check_brms(pa_sub_pressure, log_pp_check = TRUE)
  DHARMa.check_brms.all(pa_sub_pressure, integer = TRUE, outliers_type = 'bootstrap')
}
## Model has no intercept. VIFs may not be sensible.
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                      Term  VIF   VIF 95% CI Increased SE Tolerance
##                      is_A 2.33 [2.23, 2.43]         1.53      0.43
##                      is_B 1.06 [1.03, 1.10]         1.03      0.95
##     is_A:pressure_self_cw 1.02 [1.00, 1.09]         1.01      0.98
##     is_B:pressure_self_cw 1.02 [1.00, 1.09]         1.01      0.98
##  is_A:pressure_partner_cw 1.02 [1.00, 1.09]         1.01      0.98
##  is_B:pressure_partner_cw 3.70 [3.53, 3.88]         1.92      0.27
##     is_A:pressure_self_cb 4.13 [3.94, 4.34]         2.03      0.24
##  Tolerance 95% CI
##      [0.41, 0.45]
##      [0.91, 0.97]
##      [0.92, 1.00]
##      [0.92, 1.00]
##      [0.92, 1.00]
##      [0.26, 0.28]
##      [0.23, 0.25]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        cauchy         94%
##     lognormal          6%
## 
## Predicted Distribution of Response
## 
##                Distribution Probability
##  neg. binomial (zero-infl.)         78%
##               beta-binomial         16%
##                   lognormal          3%
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 13 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching
## for problematic observations.

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 6, observations = 3736, p-value = 1
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0.0005353319 0.0029443255
## sample estimates:
## outlier frequency (expected: 0.00158458244111349 ) 
##                                        0.001605996
# rope range for continuous part of the model
rope_factor <- sd(log(pa_sub_pressure$data$pa_sub[pa_sub_pressure$data$pa_sub > 0]))
rope_range_continuous = c(-0.1 * rope_factor, 0.1 * rope_factor)

summary_pa_sub <- summarize_brms(
  pa_sub_pressure, 
  stats_to_report = stats_to_report,
  rope_range = rope_range_continuous,
  hu_rope_range = c(-0.18, 0.18),
  #model_rows_fixed = model_rows_fixed,
  #model_rows_random = model_rows_random,
  #model_rownames_fixed = model_rownames_fixed,
  #model_rownames_random = model_rownames_random,
  exponentiate = T) 
## Model has no intercept. VIFs may not be sensible.
## Model has no intercept. VIFs may not be sensible.
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior
##   samples is recommended.
## Warning in summarize_brms(pa_sub_pressure, stats_to_report =
## stats_to_report, : Coefficients were exponentiated. Double check if this
## was intended.
# Print the updated dataframe
summary_pa_sub %>%
  print_df(
    #rows_to_pack = rows_to_pack
    )
exp(Est.)_hu SE_hu 95% CI_hu pd_hu ROPE_hu inside ROPE_hu BF_hu BF_Evidence_hu Rhat_hu Bulk_ESS_hu Tail_ESS_hu exp(Est.)_nonzero SE_nonzero 95% CI_nonzero pd_nonzero ROPE_nonzero inside ROPE_nonzero BF_nonzero BF_Evidence_nonzero Rhat_nonzero Bulk_ESS_nonzero Tail_ESS_nonzero
is_A 0.86 0.16 [ 0.62, 1.22] 0.778 [0.84, 1.20] 0.543 0.093 Strong Evidence for Null 1.002 1619 2259 48.30*** 3.53 [41.60, 55.61] 1.000 [0.92, 1.08] 0.000 >100 Overwhelming Evidence 1.001 1078 1836
is_A:day 0.68* 0.11 [ 0.49, 0.94] 0.992 [0.84, 1.20] 0.102 1.264 Weak Evidence 1.000 3627 2877 1.09 0.09 [ 0.93, 1.27] 0.835 [0.92, 1.08] 0.455 0.011 Very Strong Evidence for Null 1.000 5148 3217
is_A:pressure_partner_cb 0.63 0.38 [ 0.19, 2.19] 0.768 [0.84, 1.20] 0.170 0.441 Weak Evidence for Null 1.002 1387 2231 0.94 0.27 [ 0.54, 1.65] 0.581 [0.92, 1.08] 0.210 0.017 Very Strong Evidence for Null 1.004 1193 2193
is_A:pressure_partner_cw 2.43*** 0.57 [ 1.58, 4.20] 1.000 [0.84, 1.20] 0.000 >100 Overwhelming Evidence 1.001 3014 2270 0.94 0.05 [ 0.84, 1.05] 0.882 [0.92, 1.08] 0.574 0.013 Very Strong Evidence for Null 1.000 3919 2837
is_A:pressure_self_cb 0.56 0.43 [ 0.13, 2.53] 0.774 [0.84, 1.20] 0.143 0.497 Weak Evidence for Null 1.001 1507 2224 1.56 0.52 [ 0.79, 3.18] 0.907 [0.92, 1.08] 0.071 0.040 Strong Evidence for Null 1.005 1500 2130
is_A:pressure_self_cw 1.82* 0.59 [ 1.03, 4.33] 0.978 [0.84, 1.20] 0.061 1.476 Weak Evidence 1.001 2051 1761 0.85* 0.06 [ 0.72, 0.99] 0.980 [0.92, 1.08] 0.135 0.079 Strong Evidence for Null 1.001 2715 2077
is_B 1.14 0.22 [ 0.80, 1.65] 0.769 [0.84, 1.20] 0.551 0.095 Strong Evidence for Null 1.002 1294 1946 52.19*** 4.31 [44.52, 61.49] 1.000 [0.92, 1.08] 0.000 >100 Overwhelming Evidence 1.003 875 1955
is_B:day 0.56*** 0.10 [ 0.39, 0.79] 1.000 [0.84, 1.20] 0.011 17.113 Strong Evidence 1.000 4000 3081 0.83* 0.07 [ 0.72, 0.97] 0.990 [0.92, 1.08] 0.090 0.109 Moderate Evidence for Null 1.000 4854 3016
is_B:pressure_partner_cb 1.15 0.87 [ 0.26, 5.29] 0.568 [0.84, 1.20] 0.186 0.357 Weak Evidence for Null 1.001 1340 1906 1.09 0.44 [ 0.48, 2.50] 0.586 [0.92, 1.08] 0.156 0.020 Very Strong Evidence for Null 1.003 1237 1738
is_B:pressure_partner_cw 3.33* 2.03 [ 1.23, 16.18] 0.989 [0.84, 1.20] 0.017 4.570 Moderate Evidence 1.001 1929 2039 0.99 0.07 [ 0.82, 1.13] 0.567 [0.92, 1.08] 0.718 0.007 Very Strong Evidence for Null 1.000 2374 2131
is_B:pressure_self_cb 0.48 0.29 [ 0.13, 1.54] 0.893 [0.84, 1.20] 0.116 0.615 Weak Evidence for Null 1.001 1311 2034 1.08 0.34 [ 0.55, 2.03] 0.591 [0.92, 1.08] 0.189 0.019 Very Strong Evidence for Null 1.001 1278 2026
is_B:pressure_self_cw 1.29 0.35 [ 0.65, 2.65] 0.828 [0.84, 1.20] 0.316 0.229 Moderate Evidence for Null 1.001 2208 1848 0.94 0.06 [ 0.82, 1.08] 0.834 [0.92, 1.08] 0.548 0.012 Very Strong Evidence for Null 1.002 2727 2429
sd(is_A) 0.88 0.12 [0.69, 1.15] NA NA NA NA NA 1.001 1982 2636 0.32 0.05 [0.24, 0.43] NA NA NA NA NA 1.004 1243 2039
sd(is_A:pressure_partner_cw) 0.30 0.28 [0.01, 1.22] NA NA NA NA NA 1.001 1672 1800 0.06 0.05 [0.00, 0.21] NA NA NA NA NA 1.000 2120 2016
sd(is_A:pressure_self_cw) 0.70 0.51 [0.06, 2.24] NA NA NA NA NA 1.004 1258 1332 0.10 0.08 [0.00, 0.37] NA NA NA NA NA 1.001 1512 1566
sd(is_B) 0.91 0.12 [0.70, 1.19] NA NA NA NA NA 1.002 1828 2298 0.40 0.06 [0.31, 0.52] NA NA NA NA NA 1.005 1150 2250
sd(is_B:pressure_partner_cw) 1.67 0.87 [0.36, 3.91] NA NA NA NA NA 1.002 1273 960 0.09 0.08 [0.00, 0.37] NA NA NA NA NA 1.001 1324 1668
sd(is_B:pressure_self_cw) 0.75 0.49 [0.06, 2.09] NA NA NA NA NA 1.005 1011 980 0.09 0.07 [0.01, 0.28] NA NA NA NA NA 1.001 1633 1710
sigma NA NA NA NA NA NA NA NA NA NA NA 0.69 0.01 [0.67, 0.72] NA NA NA NA NA 1.000 4519 3405

Hurdle part of the model on the left, non-zero part towards the right side of the table

Pushing

formula <- bf(
  pa_sub ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:pushing_self_cw + is_B:pushing_self_cw + 
    is_A:pushing_partner_cw + is_B:pushing_partner_cw +
    
    is_A:pushing_self_cb + is_B:pushing_self_cb + 
    is_A:pushing_partner_cb + is_B:pushing_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    barriers_self_cw + barriers_self_cb + 
    is_A:day + is_B:day + 
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:pushing_self_cw + is_B:pushing_self_cw + 
    is_A:pushing_partner_cw + is_B:pushing_partner_cw |dd| coupleID),
  
  hu = ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:pushing_self_cw + is_B:pushing_self_cw + 
    is_A:pushing_partner_cw + is_B:pushing_partner_cw +
    
    is_A:pushing_self_cb + is_B:pushing_self_cb + 
    is_A:pushing_partner_cb + is_B:pushing_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    is_A:day + is_B:day + 
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:pushing_self_cw + is_B:pushing_self_cw + 
    is_A:pushing_partner_cw + is_B:pushing_partner_cw |dd| coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
) 

prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 2)", class = "b", dpar = "hu")
  
  , brms::set_prior("normal(0, 50)", class = "b", coef = 'is_A')
  , brms::set_prior("normal(0, 50)", class = "b", coef = 'is_B')
  
  , brms::set_prior("normal(0.5, 2.5)", class = "b", coef = 'is_A', dpar = 'hu')
  , brms::set_prior("normal(0.5, 2.5)", class = "b", coef = 'is_B', dpar = 'hu')
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = hurdle_lognormal()
#)

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_sub_pushing <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::hurdle_lognormal(), 
  #family = brms::hurdle_negbinomial(), 
  #family = brms::hurdle_poisson(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 42,
  file = file.path("models_cache_brms", paste0("pa_sub_hu_lognormal_pushing", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
if (check_models) {
  check_brms(pa_sub_pushing, log_pp_check = TRUE)
  DHARMa.check_brms.all(pa_sub_pushing, integer = TRUE, outliers_type = 'bootstrap')
}
## Model has no intercept. VIFs may not be sensible.
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                     Term  VIF   VIF 95% CI Increased SE Tolerance
##                     is_A 1.84 [1.76, 1.93]         1.36      0.54
##                     is_B 1.09 [1.06, 1.14]         1.04      0.92
##         barriers_self_cw 1.05 [1.02, 1.10]         1.02      0.95
##         barriers_self_cb 1.04 [1.02, 1.10]         1.02      0.96
##     is_A:pushing_self_cw 1.08 [1.05, 1.13]         1.04      0.92
##     is_B:pushing_self_cw 1.23 [1.19, 1.28]         1.11      0.81
##  is_A:pushing_partner_cw 1.31 [1.26, 1.37]         1.15      0.76
##  is_B:pushing_partner_cw 1.24 [1.20, 1.29]         1.11      0.81
##     is_A:pushing_self_cb 1.28 [1.23, 1.33]         1.13      0.78
##  Tolerance 95% CI
##      [0.52, 0.57]
##      [0.88, 0.94]
##      [0.91, 0.98]
##      [0.91, 0.98]
##      [0.88, 0.95]
##      [0.78, 0.84]
##      [0.73, 0.79]
##      [0.77, 0.84]
##      [0.75, 0.81]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        cauchy         91%
##     lognormal          9%
## 
## Predicted Distribution of Response
## 
##                Distribution Probability
##  neg. binomial (zero-infl.)         84%
##               beta-binomial          9%
##                   lognormal          6%
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 10 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching
## for problematic observations.

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 3, observations = 1776, p-value =
## 0.94
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0.000000000 0.003673986
## sample estimates:
## outlier frequency (expected: 0.00148085585585586 ) 
##                                        0.001689189
# rope range for continuous part of the model
rope_factor <- sd(log(pa_sub_pushing$data$pa_sub[pa_sub_pushing$data$pa_sub > 0]))
rope_range_continuous = c(-0.1 * rope_factor, 0.1 * rope_factor)

summary_pa_sub <- summarize_brms(
  pa_sub_pushing, 
  stats_to_report = stats_to_report,
  rope_range = rope_range_continuous,
  hu_rope_range = c(-0.18, 0.18),
  #model_rows_fixed = model_rows_fixed,
  #model_rows_random = model_rows_random,
  #model_rownames_fixed = model_rownames_fixed,
  #model_rownames_random = model_rownames_random,
  exponentiate = T) 
## Model has no intercept. VIFs may not be sensible.
## Model has no intercept. VIFs may not be sensible.
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior
##   samples is recommended.
## Warning in summarize_brms(pa_sub_pushing, stats_to_report =
## stats_to_report, : Coefficients were exponentiated. Double check if this
## was intended.
# Print the updated dataframe
summary_pa_sub %>%
  print_df(
    #rows_to_pack = rows_to_pack
    )
exp(Est.)_hu SE_hu 95% CI_hu pd_hu ROPE_hu inside ROPE_hu BF_hu BF_Evidence_hu Rhat_hu Bulk_ESS_hu Tail_ESS_hu exp(Est.)_nonzero SE_nonzero 95% CI_nonzero pd_nonzero ROPE_nonzero inside ROPE_nonzero BF_nonzero BF_Evidence_nonzero Rhat_nonzero Bulk_ESS_nonzero Tail_ESS_nonzero
barriers_self_cb NA NA NA NA NA NA NA NA NA NA NA 0.97 0.13 [ 0.75, 1.27] 0.605 [0.93, 1.08] 0.416 0.019 Very Strong Evidence for Null 1.001 1764 2409
barriers_self_cw NA NA NA NA NA NA NA NA NA NA NA 1.07* 0.03 [ 1.01, 1.13] 0.985 [0.93, 1.08] 0.651 0.093 Strong Evidence for Null 1.001 7429 3025
is_A 2.61*** 0.68 [ 1.57, 4.45] 1.000 [0.84, 1.20] 0.002 49.120 Very Strong Evidence 1.001 2799 3158 53.29*** 4.62 [44.75, 63.15] 1.000 [0.93, 1.08] 0.000 >100 Overwhelming Evidence 1.001 1898 2276
is_A:day 0.82 0.25 [ 0.46, 1.45] 0.748 [0.84, 1.20] 0.375 0.182 Moderate Evidence for Null 0.999 5470 3172 1.08 0.11 [ 0.89, 1.31] 0.795 [0.93, 1.08] 0.433 0.011 Very Strong Evidence for Null 1.000 5968 2947
is_A:pushing_partner_cb 0.91 0.89 [ 0.11, 6.51] 0.537 [0.84, 1.20] 0.146 0.463 Weak Evidence for Null 1.001 2860 3041 1.14 0.41 [ 0.55, 2.32] 0.638 [0.93, 1.08] 0.149 0.021 Very Strong Evidence for Null 1.002 1511 2291
is_A:pushing_partner_cw 1.96*** 0.44 [ 1.30, 3.48] 1.000 [0.84, 1.20] 0.010 19.295 Strong Evidence 1.000 3298 2482 0.96 0.03 [ 0.89, 1.03] 0.890 [0.93, 1.08] 0.833 0.019 Very Strong Evidence for Null 1.002 4654 2884
is_A:pushing_self_cb 0.39 0.34 [ 0.07, 2.33] 0.861 [0.84, 1.20] 0.088 0.808 Weak Evidence for Null 1.000 3309 2824 1.28 0.41 [ 0.69, 2.46] 0.777 [0.93, 1.08] 0.145 0.025 Very Strong Evidence for Null 1.001 1828 2493
is_A:pushing_self_cw 1.93 0.79 [ 0.94, 5.66] 0.964 [0.84, 1.20] 0.074 1.112 Weak Evidence 1.001 2503 2226 1.00 0.06 [ 0.87, 1.12] 0.530 [0.93, 1.08] 0.798 0.012 Very Strong Evidence for Null 1.001 2157 2444
is_B 4.00*** 1.26 [ 2.24, 7.32] 1.000 [0.84, 1.20] 0.000 >100 Overwhelming Evidence 1.000 3081 2528 52.70*** 4.59 [44.27, 62.31] 1.000 [0.93, 1.08] 0.000 >100 Overwhelming Evidence 1.001 1867 2678
is_B:day 0.64 0.20 [ 0.36, 1.16] 0.920 [0.84, 1.20] 0.186 0.390 Weak Evidence for Null 1.002 6963 2937 0.94 0.09 [ 0.79, 1.13] 0.759 [0.93, 1.08] 0.492 0.010 Very Strong Evidence for Null 1.001 8205 2886
is_B:pushing_partner_cb 0.79 0.78 [ 0.10, 5.57] 0.592 [0.84, 1.20] 0.142 0.472 Weak Evidence for Null 1.001 2790 2898 1.03 0.35 [ 0.52, 2.05] 0.533 [0.93, 1.08] 0.177 0.022 Very Strong Evidence for Null 1.002 2021 2507
is_B:pushing_partner_cw 2.32* 1.13 [ 1.04, 8.50] 0.979 [0.84, 1.20] 0.040 1.872 Weak Evidence 1.001 2151 2430 1.02 0.06 [ 0.91, 1.14] 0.660 [0.93, 1.08] 0.778 0.013 Very Strong Evidence for Null 1.000 3161 2928
is_B:pushing_self_cb 0.40 0.41 [ 0.05, 3.16] 0.807 [0.84, 1.20] 0.084 0.881 Weak Evidence for Null 1.001 3290 3081 1.47 0.59 [ 0.67, 3.28] 0.847 [0.93, 1.08] 0.095 0.037 Strong Evidence for Null 1.000 1905 2362
is_B:pushing_self_cw 1.16 0.25 [ 0.72, 1.90] 0.745 [0.84, 1.20] 0.477 0.158 Moderate Evidence for Null 1.001 2795 2632 0.94 0.03 [ 0.88, 1.01] 0.943 [0.93, 1.08] 0.708 0.028 Very Strong Evidence for Null 1.000 4564 3152
sd(is_A) 1.17 0.20 [0.86, 1.65] NA NA NA NA NA 1.001 2876 3069 0.36 0.05 [0.26, 0.50] NA NA NA NA NA 1.001 1765 2700
sd(is_A:pushing_partner_cw) 0.66 0.33 [0.10, 1.48] NA NA NA NA NA 1.000 1570 1332 0.04 0.04 [0.00, 0.14] NA NA NA NA NA 1.003 1972 2212
sd(is_A:pushing_self_cw) 1.54 0.59 [0.51, 3.12] NA NA NA NA NA 1.002 1619 1623 0.18 0.06 [0.06, 0.33] NA NA NA NA NA 1.001 1388 1087
sd(is_B) 1.38 0.23 [0.99, 1.92] NA NA NA NA NA 1.000 2820 2960 0.39 0.06 [0.29, 0.54] NA NA NA NA NA 1.001 1811 2685
sd(is_B:pushing_partner_cw) 1.74 0.84 [0.24, 3.69] NA NA NA NA NA 1.002 1076 922 0.17 0.07 [0.03, 0.31] NA NA NA NA NA 1.001 1494 1330
sd(is_B:pushing_self_cw) 0.85 0.32 [0.32, 1.74] NA NA NA NA NA 1.001 1714 1917 0.05 0.04 [0.00, 0.15] NA NA NA NA NA 1.001 1790 2025
sigma NA NA NA NA NA NA NA NA NA NA NA 0.67 0.01 [0.64, 0.70] NA NA NA NA NA 1.001 4994 2914

Hurdle part of the model on the left, non-zero part towards the right side of the table

Device Based MVPA

range(df_double$pa_obj, na.rm = T) 
## [1]   5.75 971.25
hist(df_double$pa_obj, breaks = 50)

df_double$pa_obj_log <- log(df_double$pa_obj)

hist(df_double$pa_obj_log, breaks = 50)

Persuasion

formula <- bf(
  pa_obj ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:persuasion_self_cw + is_B:persuasion_self_cw + 
    is_A:persuasion_partner_cw + is_B:persuasion_partner_cw +
    
    is_A:persuasion_self_cb + is_B:persuasion_self_cb + 
    is_A:persuasion_partner_cb + is_B:persuasion_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    is_A:day + is_B:day + 
    
    is_A:weartime_self_cw + is_B:weartime_self_cw + 
    is_A:weartime_self_cb + is_B:weartime_self_cb +
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:persuasion_self_cw + is_B:persuasion_self_cw + 
    is_A:persuasion_partner_cw + is_B:persuasion_partner_cw |dd| coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
)



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  
  , brms::set_prior("normal(0, 50)", class = "b", coef = 'is_A')
  , brms::set_prior("normal(0, 50)", class = "b", coef = 'is_B')
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = lognormal()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_obj_log_persuasion <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = lognormal(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("pa_obj_log_gaussian_persuasion", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
if (check_models) {
  check_brms(pa_obj_log_persuasion, log_pp_check = TRUE)
  DHARMa.check_brms.all(pa_obj_log_persuasion, integer = TRUE, outliers_type = 'bootstrap')
}
## Model has no intercept. VIFs may not be sensible.
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                        Term  VIF   VIF 95% CI Increased SE Tolerance
##                        is_A 1.94 [1.86, 2.04]         1.39      0.51
##                        is_B 1.09 [1.06, 1.14]         1.04      0.92
##     is_A:persuasion_self_cw 1.17 [1.13, 1.21]         1.08      0.86
##     is_B:persuasion_self_cw 1.12 [1.08, 1.16]         1.06      0.90
##  is_A:persuasion_partner_cw 1.05 [1.02, 1.10]         1.02      0.95
##  is_B:persuasion_partner_cw 3.62 [3.43, 3.82]         1.90      0.28
##     is_A:persuasion_self_cb 2.84 [2.70, 2.99]         1.68      0.35
##     is_B:persuasion_self_cb 3.82 [3.62, 4.04]         1.95      0.26
##  is_A:persuasion_partner_cb 2.87 [2.72, 3.02]         1.69      0.35
##  Tolerance 95% CI
##      [0.49, 0.54]
##      [0.88, 0.94]
##      [0.82, 0.88]
##      [0.86, 0.92]
##      [0.91, 0.98]
##      [0.26, 0.29]
##      [0.33, 0.37]
##      [0.25, 0.28]
##      [0.33, 0.37]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        cauchy         88%
##     lognormal          9%
##       weibull          3%
## 
## Predicted Distribution of Response
## 
##                Distribution Probability
##                   lognormal         72%
##                     tweedie         16%
##  neg. binomial (zero-infl.)          9%
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## DHARMa:testOutliers with type = binomial may have inflated Type I error rates for integer-valued distributions. To get a more exact result, it is recommended to re-run testOutliers with type = 'bootstrap'. See ?testOutliers for details

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 31, observations = 3343, p-value <
## 2.2e-16
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0.001495663 0.004786120
## sample estimates:
## outlier frequency (expected: 0.00320370924319474 ) 
##                                        0.009273108
# rope range for lognormal model
rope_factor <- sd(log(pa_obj_log_persuasion$data$pa_obj))
rope_range_log = c(-0.1 * rope_factor, 0.1 * rope_factor)

summary_pa_obj <- summarize_brms(
  pa_obj_log_persuasion, 
  stats_to_report = stats_to_report,
  rope_range = rope_range_log,
  #model_rows_fixed = model_rows_fixed,
  #model_rows_random = model_rows_random,
  #model_rownames_fixed = model_rownames_fixed,
  #model_rownames_random = model_rownames_random,
  exponentiate = T) 
## Model has no intercept. VIFs may not be sensible.
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior
##   samples is recommended.
## Warning in summarize_brms(pa_obj_log_persuasion, stats_to_report =
## stats_to_report, : Coefficients were exponentiated. Double check if this
## was intended.
summary_pa_obj %>%
  print_df(
    #rows_to_pack = rows_to_pack
    )
exp(Est.) SE 95% CI pd ROPE inside ROPE BF BF_Evidence Rhat Bulk_ESS Tail_ESS
is_A 101.71*** 6.22 [ 89.97, 115.11] 1.000 [0.94, 1.07] 0.000 >100 Overwhelming Evidence 1.003 951 1434
is_B 136.08*** 7.27 [122.13, 151.56] 1.000 [0.94, 1.07] 0.000 >100 Overwhelming Evidence 1.005 1389 1918
is_A:persuasion_self_cw 1.05* 0.02 [ 1.01, 1.09] 0.987 [0.94, 1.07] 0.843 0.076 Strong Evidence for Null 1.000 3928 3337
is_B:persuasion_self_cw 1.01 0.02 [ 0.97, 1.05] 0.655 [0.94, 1.07] 0.997 0.006 Very Strong Evidence for Null 1.001 3630 3130
is_A:persuasion_partner_cw 1.02 0.02 [ 0.98, 1.06] 0.796 [0.94, 1.07] 0.987 0.008 Very Strong Evidence for Null 0.999 3459 3176
is_B:persuasion_partner_cw 1.03 0.02 [ 0.99, 1.07] 0.913 [0.94, 1.07] 0.964 0.014 Very Strong Evidence for Null 1.000 3753 3537
is_A:persuasion_self_cb 1.25 0.22 [ 0.87, 1.77] 0.899 [0.94, 1.07] 0.128 0.036 Strong Evidence for Null 1.004 994 1743
is_B:persuasion_self_cb 0.89 0.11 [ 0.69, 1.13] 0.837 [0.94, 1.07] 0.264 0.018 Very Strong Evidence for Null 1.003 1266 1932
is_A:persuasion_partner_cb 0.89 0.14 [ 0.67, 1.21] 0.785 [0.94, 1.07] 0.244 0.021 Very Strong Evidence for Null 1.004 1106 2003
is_B:persuasion_partner_cb 1.29 0.19 [ 0.97, 1.76] 0.960 [0.94, 1.07] 0.086 0.059 Strong Evidence for Null 1.003 1153 2022
is_A:day 0.98 0.05 [ 0.90, 1.08] 0.626 [0.94, 1.07] 0.805 0.004 Very Strong Evidence for Null 1.001 9420 2738
is_B:day 0.95 0.05 [ 0.86, 1.04] 0.856 [0.94, 1.07] 0.608 0.006 Very Strong Evidence for Null 1.003 8586 2175
is_A:weartime_self_cw 1.00*** 0.00 [ 1.00, 1.00] 1.000 [0.94, 1.07] 1.000 3.072 Moderate Evidence 1.000 9905 3242
is_B:weartime_self_cw 1.00* 0.00 [ 1.00, 1.00] 0.986 [0.94, 1.07] 1.000 0.042 Strong Evidence for Null 1.002 7670 2893
is_A:weartime_self_cb 1.00 0.00 [ 1.00, 1.00] 0.702 [0.94, 1.07] 1.000 0.013 Very Strong Evidence for Null 1.002 1124 1965
is_B:weartime_self_cb 1.00 0.00 [ 1.00, 1.00] 0.575 [0.94, 1.07] 1.000 0.011 Very Strong Evidence for Null 1.001 1304 1992
sd(is_A) 0.34 0.04 [0.27, 0.44] NA NA NA NA NA 1.001 1313 2072
sd(is_B) 0.29 0.04 [0.22, 0.38] NA NA NA NA NA 1.001 1472 2002
sd(is_A:persuasion_self_cw) 0.06 0.02 [0.01, 0.11] NA NA NA NA NA 1.006 851 503
sd(is_B:persuasion_self_cw) 0.06 0.02 [0.02, 0.11] NA NA NA NA NA 1.002 1911 1806
sd(is_A:persuasion_partner_cw) 0.07 0.02 [0.02, 0.13] NA NA NA NA NA 1.001 1384 1295
sd(is_B:persuasion_partner_cw) 0.07 0.03 [0.01, 0.13] NA NA NA NA NA 1.002 1208 1198
sigma 0.55 0.01 [0.53, 0.56] NA NA NA NA NA 1.000 8451 3220

Pressure

formula <- bf(
  pa_obj ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:pressure_self_cw + is_B:pressure_self_cw + 
    is_A:pressure_partner_cw + is_B:pressure_partner_cw +
    
    is_A:pressure_self_cb + is_B:pressure_self_cb + 
    is_A:pressure_partner_cb + is_B:pressure_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    is_A:day + is_B:day + 
    
    is_A:weartime_self_cw + is_B:weartime_self_cw + 
    is_A:weartime_self_cb + is_B:weartime_self_cb +
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:pressure_self_cw + is_B:pressure_self_cw + 
    is_A:pressure_partner_cw + is_B:pressure_partner_cw |dd| coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
)



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  
  , brms::set_prior("normal(0, 50)", class = "b", coef = 'is_A')
  , brms::set_prior("normal(0, 50)", class = "b", coef = 'is_B')
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = lognormal()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_obj_log_pressure <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = lognormal(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("pa_obj_log_gaussian_pressure", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
if (check_models) {
  check_brms(pa_obj_log_pressure, log_pp_check = TRUE)
  DHARMa.check_brms.all(pa_obj_log_pressure, integer = TRUE, outliers_type = 'bootstrap')
}
## Model has no intercept. VIFs may not be sensible.
## Model has interaction terms. VIFs might be inflated.
##   You may check multicollinearity among predictors of a model without
##   interaction terms.
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                      Term  VIF     VIF 95% CI Increased SE Tolerance
##                      is_A 1.90 [ 1.82,  1.99]         1.38      0.53
##                      is_B 1.03 [ 1.01,  1.09]         1.01      0.97
##     is_A:pressure_self_cw 1.02 [ 1.00,  1.11]         1.01      0.98
##     is_B:pressure_self_cw 1.03 [ 1.01,  1.09]         1.01      0.97
##  is_A:pressure_partner_cw 1.02 [ 1.01,  1.10]         1.01      0.98
##  Tolerance 95% CI
##      [0.50, 0.55]
##      [0.92, 0.99]
##      [0.90, 1.00]
##      [0.92, 0.99]
##      [0.91, 0.99]
## 
## Moderate Correlation
## 
##                      Term  VIF     VIF 95% CI Increased SE Tolerance
##     is_A:pressure_self_cb 8.76 [ 8.26,  9.29]         2.96      0.11
##  is_A:pressure_partner_cb 9.02 [ 8.51,  9.57]         3.00      0.11
##  Tolerance 95% CI
##      [0.11, 0.12]
##      [0.10, 0.12]
## 
## High Correlation
## 
##                      Term   VIF     VIF 95% CI Increased SE Tolerance
##  is_B:pressure_partner_cw 10.59 [ 9.98, 11.24]         3.25      0.09
##     is_B:pressure_self_cb 10.72 [10.11, 11.38]         3.27      0.09
##  Tolerance 95% CI
##      [0.09, 0.10]
##      [0.09, 0.10]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        cauchy         88%
##     lognormal          9%
##       weibull          3%
## 
## Predicted Distribution of Response
## 
##                Distribution Probability
##                   lognormal         72%
##                     tweedie         16%
##  neg. binomial (zero-infl.)          9%
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 3 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching
## for problematic observations.
## DHARMa:testOutliers with type = binomial may have inflated Type I error rates for integer-valued distributions. To get a more exact result, it is recommended to re-run testOutliers with type = 'bootstrap'. See ?testOutliers for details

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 29, observations = 3337, p-value <
## 2.2e-16
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0.001498352 0.005094396
## sample estimates:
## outlier frequency (expected: 0.00310758166017381 ) 
##                                        0.008690441
# rope range for lognormal model
rope_factor <- sd(log(pa_obj_log_pressure$data$pa_obj))
rope_range_log = c(-0.1 * rope_factor, 0.1 * rope_factor)

summary_pa_obj <- summarize_brms(
  pa_obj_log_pressure, 
  stats_to_report = stats_to_report,
  rope_range = rope_range_log,
  #model_rows_fixed = model_rows_fixed,
  #model_rows_random = model_rows_random,
  #model_rownames_fixed = model_rownames_fixed,
  #model_rownames_random = model_rownames_random,
  exponentiate = T) 
## Model has no intercept. VIFs may not be sensible.
## Model has interaction terms. VIFs might be inflated.
##   You may check multicollinearity among predictors of a model without
##   interaction terms.
## Warning in compute_rope(range = rope_range): Collinearity detected. Some
## VIFs are > 10. This may invalidate ROPE inferences!

Check for Multicollinearity

Low Correlation

                 Term  VIF     VIF 95% CI Increased SE Tolerance
                 is_A 1.90 [ 1.82,  1.99]         1.38      0.53
                 is_B 1.03 [ 1.01,  1.09]         1.01      0.97
is_A:pressure_self_cw 1.02 [ 1.00,  1.11]         1.01      0.98
is_B:pressure_self_cw 1.03 [ 1.01,  1.09]         1.01      0.97

is_A:pressure_partner_cw 1.02 [ 1.01, 1.10] 1.01 0.98 Tolerance 95% CI [0.50, 0.55] [0.92, 0.99] [0.90, 1.00] [0.92, 0.99] [0.91, 0.99]

Moderate Correlation

                 Term  VIF     VIF 95% CI Increased SE Tolerance
is_A:pressure_self_cb 8.76 [ 8.26,  9.29]         2.96      0.11

is_A:pressure_partner_cb 9.02 [ 8.51, 9.57] 3.00 0.11 Tolerance 95% CI [0.11, 0.12] [0.10, 0.12]

High Correlation

                 Term   VIF     VIF 95% CI Increased SE Tolerance

is_B:pressure_partner_cw 10.59 [ 9.98, 11.24] 3.25 0.09 is_B:pressure_self_cb 10.72 [10.11, 11.38] 3.27 0.09 Tolerance 95% CI [0.09, 0.10] [0.09, 0.10]

## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior
##   samples is recommended.
## Warning in summarize_brms(pa_obj_log_pressure, stats_to_report =
## stats_to_report, : Coefficients were exponentiated. Double check if this
## was intended.
summary_pa_obj %>%
  print_df(
    #rows_to_pack = rows_to_pack
    )
exp(Est.) SE 95% CI pd ROPE inside ROPE BF BF_Evidence Rhat Bulk_ESS Tail_ESS
is_A 101.82*** 6.29 [ 90.01, 115.55] 1.000 [0.94, 1.07] 0.000 >100 Overwhelming Evidence 1.001 760 1414
is_B 137.46*** 7.80 [123.29, 153.89] 1.000 [0.94, 1.07] 0.000 >100 Overwhelming Evidence 1.001 939 2113
is_A:pressure_self_cw 0.99 0.05 [ 0.87, 1.09] 0.607 [0.94, 1.07] 0.760 0.005 Very Strong Evidence for Null 1.001 3475 2940
is_B:pressure_self_cw 0.99 0.04 [ 0.91, 1.08] 0.609 [0.94, 1.07] 0.853 0.005 Very Strong Evidence for Null 0.999 5396 3072
is_A:pressure_partner_cw 1.01 0.04 [ 0.93, 1.09] 0.600 [0.94, 1.07] 0.898 0.004 Very Strong Evidence for Null 1.000 5928 3035
is_B:pressure_partner_cw 1.02 0.06 [ 0.92, 1.16] 0.665 [0.94, 1.07] 0.713 0.006 Very Strong Evidence for Null 1.001 3788 2687
is_A:pressure_self_cb 0.85 0.29 [ 0.43, 1.66] 0.690 [0.94, 1.07] 0.132 0.017 Very Strong Evidence for Null 1.004 1534 1963
is_B:pressure_self_cb 0.93 0.23 [ 0.58, 1.53] 0.618 [0.94, 1.07] 0.202 0.012 Very Strong Evidence for Null 1.003 1593 2174
is_A:pressure_partner_cb 1.24 0.34 [ 0.73, 2.15] 0.791 [0.94, 1.07] 0.141 0.021 Very Strong Evidence for Null 1.005 1437 2207
is_B:pressure_partner_cb 1.21 0.38 [ 0.64, 2.23] 0.726 [0.94, 1.07] 0.123 0.016 Very Strong Evidence for Null 1.004 1678 2434
is_A:day 0.95 0.04 [ 0.87, 1.04] 0.851 [0.94, 1.07] 0.642 0.006 Very Strong Evidence for Null 1.002 8015 3016
is_B:day 0.92 0.04 [ 0.84, 1.01] 0.964 [0.94, 1.07] 0.352 0.019 Very Strong Evidence for Null 1.003 9930 2577
is_A:weartime_self_cw 1.00*** 0.00 [ 1.00, 1.00] 1.000 [0.94, 1.07] 1.000 3.889 Moderate Evidence 1.000 7074 2735
is_B:weartime_self_cw 1.00 0.00 [ 1.00, 1.00] 0.971 [0.94, 1.07] 1.000 0.026 Very Strong Evidence for Null 1.001 8048 2482
is_A:weartime_self_cb 1.00 0.00 [ 1.00, 1.00] 0.626 [0.94, 1.07] 1.000 0.011 Very Strong Evidence for Null 1.003 914 1459
is_B:weartime_self_cb 1.00 0.00 [ 1.00, 1.00] 0.718 [0.94, 1.07] 1.000 0.012 Very Strong Evidence for Null 1.009 1187 2234
sd(is_A) 0.34 0.04 [0.27, 0.44] NA NA NA NA NA 1.003 1060 1993
sd(is_B) 0.29 0.04 [0.23, 0.39] NA NA NA NA NA 1.002 1455 2240
sd(is_A:pressure_self_cw) 0.08 0.07 [0.00, 0.25] NA NA NA NA NA 1.000 1310 1811
sd(is_B:pressure_self_cw) 0.05 0.04 [0.00, 0.18] NA NA NA NA NA 1.001 1935 1832
sd(is_A:pressure_partner_cw) 0.04 0.03 [0.00, 0.13] NA NA NA NA NA 1.002 2463 1819
sd(is_B:pressure_partner_cw) 0.07 0.06 [0.00, 0.26] NA NA NA NA NA 1.000 1722 1765
sigma 0.56 0.01 [0.54, 0.57] NA NA NA NA NA 1.001 7882 2946

Pushing

formula <- bf(
  pa_obj ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:pushing_self_cw + is_B:pushing_self_cw + 
    is_A:pushing_partner_cw + is_B:pushing_partner_cw +
    
    is_A:pushing_self_cb + is_B:pushing_self_cb + 
    is_A:pushing_partner_cb + is_B:pushing_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    barriers_self_cw + barriers_self_cb + 
    is_A:day + is_B:day + 
    
    is_A:weartime_self_cw + is_B:weartime_self_cw + 
    is_A:weartime_self_cb + is_B:weartime_self_cb +
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:pushing_self_cw + is_B:pushing_self_cw + 
    is_A:pushing_partner_cw + is_B:pushing_partner_cw |dd| coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
)



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  
  , brms::set_prior("normal(0, 50)", class = "b", coef = 'is_A')
  , brms::set_prior("normal(0, 50)", class = "b", coef = 'is_B')
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = lognormal()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_obj_log_pushing <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = lognormal(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("pa_obj_log_gaussian_pushing", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
if (check_models) {
  check_brms(pa_obj_log_pushing, log_pp_check = TRUE)
  DHARMa.check_brms.all(pa_obj_log_pushing, integer = TRUE, outliers_type = 'bootstrap')
}
## Model has no intercept. VIFs may not be sensible.
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                     Term  VIF   VIF 95% CI Increased SE Tolerance
##                     is_A 2.03 [1.92, 2.17]         1.43      0.49
##                     is_B 1.02 [1.00, 1.14]         1.01      0.98
##         barriers_self_cw 1.09 [1.05, 1.15]         1.04      0.92
##         barriers_self_cb 1.07 [1.03, 1.13]         1.03      0.94
##     is_A:pushing_self_cw 1.06 [1.03, 1.13]         1.03      0.94
##     is_B:pushing_self_cw 1.05 [1.02, 1.12]         1.02      0.96
##  is_A:pushing_partner_cw 1.10 [1.06, 1.16]         1.05      0.91
##  is_B:pushing_partner_cw 1.91 [1.80, 2.03]         1.38      0.52
##     is_A:pushing_self_cb 1.82 [1.72, 1.94]         1.35      0.55
##     is_B:pushing_self_cb 1.93 [1.82, 2.05]         1.39      0.52
##  is_A:pushing_partner_cb 1.80 [1.70, 1.91]         1.34      0.56
##  Tolerance 95% CI
##      [0.46, 0.52]
##      [0.88, 1.00]
##      [0.87, 0.95]
##      [0.88, 0.97]
##      [0.89, 0.97]
##      [0.89, 0.98]
##      [0.86, 0.94]
##      [0.49, 0.56]
##      [0.52, 0.58]
##      [0.49, 0.55]
##      [0.52, 0.59]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        cauchy         88%
##     lognormal          9%
##       weibull          3%
## 
## Predicted Distribution of Response
## 
##                Distribution Probability
##                   lognormal         72%
##                     tweedie         16%
##  neg. binomial (zero-infl.)          9%
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 3 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching
## for problematic observations.
## DHARMa:testOutliers with type = binomial may have inflated Type I error rates for integer-valued distributions. To get a more exact result, it is recommended to re-run testOutliers with type = 'bootstrap'. See ?testOutliers for details

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 14, observations = 1594, p-value <
## 2.2e-16
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0.0006273526 0.0062735257
## sample estimates:
## outlier frequency (expected: 0.00330614805520703 ) 
##                                        0.008782936
# rope range for lognormal model
rope_factor <- sd(log(pa_obj_log_pushing$data$pa_obj))
rope_range_log = c(-0.1 * rope_factor, 0.1 * rope_factor)

summary_pa_obj <- summarize_brms(
  pa_obj_log_pushing, 
  stats_to_report = stats_to_report,
  rope_range = rope_range_log,
  #model_rows_fixed = model_rows_fixed,
  #model_rows_random = model_rows_random,
  #model_rownames_fixed = model_rownames_fixed,
  #model_rownames_random = model_rownames_random,
  exponentiate = T) 
## Model has no intercept. VIFs may not be sensible.
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior
##   samples is recommended.
## Warning in summarize_brms(pa_obj_log_pushing, stats_to_report =
## stats_to_report, : Coefficients were exponentiated. Double check if this
## was intended.
summary_pa_obj %>%
  print_df(
    #rows_to_pack = rows_to_pack
    )
exp(Est.) SE 95% CI pd ROPE inside ROPE BF BF_Evidence Rhat Bulk_ESS Tail_ESS
is_A 106.16*** 7.66 [ 93.32, 122.54] 1.000 [0.94, 1.06] 0.000 >100 Overwhelming Evidence 1.003 1682 2179
is_B 146.96*** 10.28 [127.91, 167.84] 1.000 [0.94, 1.06] 0.000 >100 Overwhelming Evidence 1.001 2191 2869
barriers_self_cw 0.97 0.02 [ 0.93, 1.00] 0.970 [0.94, 1.06] 0.923 0.032 Strong Evidence for Null 1.001 9179 3053
barriers_self_cb 0.86 0.08 [ 0.72, 1.05] 0.928 [0.94, 1.06] 0.185 0.045 Strong Evidence for Null 1.001 2102 2746
is_A:pushing_self_cw 1.07 0.04 [ 0.99, 1.15] 0.958 [0.94, 1.06] 0.460 0.043 Strong Evidence for Null 1.001 4454 3162
is_B:pushing_self_cw 1.00 0.03 [ 0.94, 1.07] 0.508 [0.94, 1.06] 0.949 0.008 Very Strong Evidence for Null 1.000 4795 3193
is_A:pushing_partner_cw 0.98 0.02 [ 0.94, 1.04] 0.744 [0.94, 1.06] 0.963 0.007 Very Strong Evidence for Null 1.000 7078 2971
is_B:pushing_partner_cw 1.06 0.04 [ 0.99, 1.14] 0.949 [0.94, 1.06] 0.513 0.030 Strong Evidence for Null 1.000 5136 3496
is_A:pushing_self_cb 0.99 0.30 [ 0.55, 1.87] 0.511 [0.94, 1.06] 0.158 0.016 Very Strong Evidence for Null 1.002 1901 2345
is_B:pushing_self_cb 0.97 0.30 [ 0.53, 1.85] 0.538 [0.94, 1.06] 0.152 0.017 Very Strong Evidence for Null 1.002 1552 2479
is_A:pushing_partner_cb 1.19 0.38 [ 0.62, 2.26] 0.699 [0.94, 1.06] 0.136 0.019 Very Strong Evidence for Null 1.001 1502 2456
is_B:pushing_partner_cb 1.42 0.43 [ 0.79, 2.60] 0.875 [0.94, 1.06] 0.090 0.032 Strong Evidence for Null 1.002 2238 2707
is_A:day 1.01 0.06 [ 0.88, 1.15] 0.527 [0.94, 1.06] 0.659 0.005 Very Strong Evidence for Null 1.000 9192 3174
is_B:day 0.89 0.06 [ 0.78, 1.01] 0.960 [0.94, 1.06] 0.190 0.030 Very Strong Evidence for Null 1.004 9607 2369
is_A:weartime_self_cw 1.00** 0.00 [ 1.00, 1.00] 0.997 [0.94, 1.06] 1.000 0.254 Moderate Evidence for Null 1.000 7285 3602
is_B:weartime_self_cw 1.00* 0.00 [ 1.00, 1.00] 0.991 [0.94, 1.06] 1.000 0.102 Moderate Evidence for Null 1.000 8180 2831
is_A:weartime_self_cb 1.00 0.00 [ 1.00, 1.00] 0.677 [0.94, 1.06] 1.000 0.014 Very Strong Evidence for Null 1.003 1643 2402
is_B:weartime_self_cb 1.00 0.00 [ 1.00, 1.00] 0.875 [0.94, 1.06] 1.000 0.028 Very Strong Evidence for Null 1.000 2133 2900
sd(is_A) 0.34 0.05 [0.26, 0.46] NA NA NA NA NA 1.004 1819 2307
sd(is_B) 0.34 0.05 [0.26, 0.45] NA NA NA NA NA 1.001 1875 2972
sd(is_A:pushing_self_cw) 0.10 0.05 [0.01, 0.22] NA NA NA NA NA 1.002 1096 1604
sd(is_B:pushing_self_cw) 0.09 0.04 [0.02, 0.17] NA NA NA NA NA 1.001 1308 1555
sd(is_A:pushing_partner_cw) 0.03 0.03 [0.00, 0.10] NA NA NA NA NA 1.000 2675 2695
sd(is_B:pushing_partner_cw) 0.10 0.04 [0.01, 0.18] NA NA NA NA NA 1.001 1561 1168
sigma 0.52 0.01 [0.50, 0.54] NA NA NA NA NA 1.001 6494 2778

Affect

range(df_double$aff, na.rm = T)
## [1] 0 5
hist(df_double$aff, breaks = 15)

Persuasion

formula <- bf(
  aff ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:persuasion_self_cw + is_B:persuasion_self_cw + 
    is_A:persuasion_partner_cw + is_B:persuasion_partner_cw +
    
    is_A:persuasion_self_cb + is_B:persuasion_self_cb + 
    is_A:persuasion_partner_cb + is_B:persuasion_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    is_A:day + is_B:day + 
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:pushing_self_cw + is_B:pushing_self_cw + 
    is_A:pushing_partner_cw + is_B:pushing_partner_cw |dd| coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
)


prior1 <- c(
  brms::set_prior("normal(0, 5)", class = "b")
  
  , brms::set_prior("normal(0, 20)", class = "b", coef = 'is_A')
  , brms::set_prior("normal(0, 20)", class = "b", coef = 'is_B')
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = gaussian()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

mood_gauss_persuasion <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = gaussian(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("mood_gauss_persuasion", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
if (check_models) {
  check_brms(mood_gauss_persuasion, log_pp_check = FALSE)
  DHARMa.check_brms.all(mood_gauss_persuasion, integer = FALSE)
}
## Model has no intercept. VIFs may not be sensible.
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                        Term  VIF   VIF 95% CI Increased SE Tolerance
##                        is_A 1.15 [1.11, 1.19]         1.07      0.87
##                        is_B 1.04 [1.02, 1.09]         1.02      0.96
##     is_A:persuasion_self_cw 1.06 [1.04, 1.11]         1.03      0.94
##     is_B:persuasion_self_cw 1.10 [1.07, 1.14]         1.05      0.91
##  is_A:persuasion_partner_cw 1.04 [1.02, 1.09]         1.02      0.96
##  is_B:persuasion_partner_cw 2.66 [2.54, 2.80]         1.63      0.38
##     is_A:persuasion_self_cb 2.69 [2.57, 2.83]         1.64      0.37
##  Tolerance 95% CI
##      [0.84, 0.90]
##      [0.92, 0.98]
##      [0.90, 0.96]
##      [0.88, 0.94]
##      [0.92, 0.98]
##      [0.36, 0.39]
##      [0.35, 0.39]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        cauchy         47%
##        normal         34%
##   exponential          6%
## 
## Predicted Distribution of Response
## 
##   Distribution Probability
##        tweedie         41%
##  beta-binomial         38%
##    half-cauchy          6%
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## 
##  DHARMa outlier test based on exact binomial test with approximate
##  expectations
## 
## data:  model.check
## outliers at both margin(s) = 38, observations = 3688, p-value =
## 1.198e-15
## alternative hypothesis: true probability of success is not equal to 0.001998002
## 95 percent confidence interval:
##  0.007301532 0.014115438
## sample estimates:
## frequency of outliers (expected: 0.001998001998002 ) 
##                                           0.01030369
summary_mood <- summarize_brms(
  mood_gauss_persuasion, 
  stats_to_report = stats_to_report,
  #model_rows_fixed = model_rows_fixed,
  #model_rows_random = model_rows_random,
  #model_rownames_fixed = model_rownames_fixed,
  #model_rownames_random = model_rownames_random,
  exponentiate = F) 
## Model has no intercept. VIFs may not be sensible.
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior
##   samples is recommended.
summary_mood %>%
  print_df(
    #rows_to_pack = rows_to_pack
    )
Est. SE 95% CI pd ROPE inside ROPE BF BF_Evidence Rhat Bulk_ESS Tail_ESS
is_A 3.61*** 0.12 [ 3.37, 3.84] 1.000 [-0.11, 0.11] 0.000 >100 Overwhelming Evidence 1.002 453 899
is_B 3.87*** 0.14 [ 3.61, 4.16] 1.000 [-0.11, 0.11] 0.000 >100 Overwhelming Evidence 1.010 446 1023
is_A:persuasion_self_cw 0.01 0.02 [-0.04, 0.05] 0.634 [-0.11, 0.11] 1.000 0.003 Very Strong Evidence for Null 1.000 5822 3143
is_B:persuasion_self_cw -0.01 0.02 [-0.06, 0.03] 0.708 [-0.11, 0.11] 1.000 0.003 Very Strong Evidence for Null 1.001 5420 2727
is_A:persuasion_partner_cw 0.05* 0.02 [ 0.00, 0.09] 0.980 [-0.11, 0.11] 0.998 0.025 Very Strong Evidence for Null 1.000 5603 2712
is_B:persuasion_partner_cw 0.01 0.02 [-0.04, 0.05] 0.633 [-0.11, 0.11] 1.000 0.003 Very Strong Evidence for Null 1.000 5288 3058
is_A:persuasion_self_cb -0.27 0.32 [-0.90, 0.39] 0.797 [-0.11, 0.11] 0.190 0.023 Very Strong Evidence for Null 1.002 571 957
is_B:persuasion_self_cb 0.38 0.33 [-0.31, 1.02] 0.868 [-0.11, 0.11] 0.142 0.036 Strong Evidence for Null 1.003 672 1116
is_A:persuasion_partner_cb 0.39 0.27 [-0.18, 0.92] 0.917 [-0.11, 0.11] 0.120 0.040 Strong Evidence for Null 1.003 583 1141
is_B:persuasion_partner_cb -0.20 0.39 [-0.97, 0.58] 0.701 [-0.11, 0.11] 0.206 0.022 Very Strong Evidence for Null 1.005 672 1191
is_A:day 0.45*** 0.07 [ 0.31, 0.59] 1.000 [-0.11, 0.11] 0.000 >100 Overwhelming Evidence 1.000 4830 2774
is_B:day 0.08 0.07 [-0.05, 0.22] 0.874 [-0.11, 0.11] 0.656 0.005 Very Strong Evidence for Null 1.000 5104 3328
sd(is_A) 0.70 0.09 [0.55, 0.90] NA NA NA NA NA 1.011 605 1204
sd(is_B) 0.79 0.10 [0.63, 1.03] NA NA NA NA NA 1.009 774 1225
sd(is_A:pushing_self_cw) 0.11 0.06 [0.01, 0.24] NA NA NA NA NA 1.009 1039 1100
sd(is_B:pushing_self_cw) 0.06 0.05 [0.00, 0.18] NA NA NA NA NA 1.002 1271 1356
sd(is_A:pushing_partner_cw) 0.13 0.05 [0.03, 0.24] NA NA NA NA NA 1.007 732 301
sd(is_B:pushing_partner_cw) 0.12 0.06 [0.01, 0.25] NA NA NA NA NA 1.003 1411 1092
sigma 0.87 0.01 [0.85, 0.89] NA NA NA NA NA 1.000 5092 2948

Pressure

formula <- bf(
  aff ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:pressure_self_cw + is_B:pressure_self_cw + 
    is_A:pressure_partner_cw + is_B:pressure_partner_cw +
    
    is_A:pressure_self_cb + is_B:pressure_self_cb + 
    is_A:pressure_partner_cb + is_B:pressure_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    is_A:day + is_B:day + 
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:pushing_self_cw + is_B:pushing_self_cw + 
    is_A:pushing_partner_cw + is_B:pushing_partner_cw |dd| coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
)


prior1 <- c(
  brms::set_prior("normal(0, 5)", class = "b")
  
  , brms::set_prior("normal(0, 20)", class = "b", coef = 'is_A')
  , brms::set_prior("normal(0, 20)", class = "b", coef = 'is_B')
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = gaussian()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

mood_gauss_pressure <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = gaussian(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("mood_gauss_pressure", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
if (check_models) {
  check_brms(mood_gauss_pressure, log_pp_check = FALSE)
  DHARMa.check_brms.all(mood_gauss_pressure, integer = FALSE)
}
## Model has no intercept. VIFs may not be sensible.
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                      Term  VIF   VIF 95% CI Increased SE Tolerance
##                      is_A 1.20 [1.17, 1.25]         1.10      0.83
##                      is_B 1.01 [1.00, 1.12]         1.01      0.99
##     is_A:pressure_self_cw 1.02 [1.00, 1.10]         1.01      0.98
##     is_B:pressure_self_cw 1.03 [1.01, 1.09]         1.01      0.97
##  is_A:pressure_partner_cw 1.00 [1.00, 3.28]         1.00      1.00
##  Tolerance 95% CI
##      [0.80, 0.86]
##      [0.89, 1.00]
##      [0.91, 1.00]
##      [0.92, 0.99]
##      [0.30, 1.00]
## 
## Moderate Correlation
## 
##                      Term  VIF   VIF 95% CI Increased SE Tolerance
##  is_B:pressure_partner_cw 8.85 [8.37, 9.37]         2.98      0.11
##     is_A:pressure_self_cb 8.07 [7.63, 8.53]         2.84      0.12
##  Tolerance 95% CI
##      [0.11, 0.12]
##      [0.12, 0.13]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        cauchy         47%
##        normal         34%
##   exponential          6%
## 
## Predicted Distribution of Response
## 
##   Distribution Probability
##        tweedie         41%
##  beta-binomial         38%
##    half-cauchy          6%
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## 
##  DHARMa outlier test based on exact binomial test with approximate
##  expectations
## 
## data:  model.check
## outliers at both margin(s) = 38, observations = 3684, p-value =
## 1.158e-15
## alternative hypothesis: true probability of success is not equal to 0.001998002
## 95 percent confidence interval:
##  0.007309471 0.014130735
## sample estimates:
## frequency of outliers (expected: 0.001998001998002 ) 
##                                           0.01031488
summary_mood <- summarize_brms(
  mood_gauss_pressure, 
  stats_to_report = stats_to_report,
  #model_rows_fixed = model_rows_fixed,
  #model_rows_random = model_rows_random,
  #model_rownames_fixed = model_rownames_fixed,
  #model_rownames_random = model_rownames_random,
  exponentiate = F) 
## Model has no intercept. VIFs may not be sensible.
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior
##   samples is recommended.
summary_mood %>%
  print_df(
    #rows_to_pack = rows_to_pack
    )
Est. SE 95% CI pd ROPE inside ROPE BF BF_Evidence Rhat Bulk_ESS Tail_ESS
is_A 3.63*** 0.14 [ 3.35, 3.90] 1.000 [-0.11, 0.11] 0.000 >100 Overwhelming Evidence 1.020 198 409
is_B 3.87*** 0.14 [ 3.58, 4.16] 1.000 [-0.11, 0.11] 0.000 >100 Overwhelming Evidence 1.018 233 687
is_A:pressure_self_cw -0.06 0.06 [-0.17, 0.06] 0.844 [-0.11, 0.11] 0.841 0.005 Very Strong Evidence for Null 1.001 6955 3236
is_B:pressure_self_cw 0.01 0.06 [-0.11, 0.11] 0.542 [-0.11, 0.11] 0.960 0.003 Very Strong Evidence for Null 1.001 4185 3236
is_A:pressure_partner_cw 0.05 0.06 [-0.06, 0.16] 0.809 [-0.11, 0.11] 0.877 0.004 Very Strong Evidence for Null 1.002 3919 3108
is_B:pressure_partner_cw 0.04 0.06 [-0.07, 0.15] 0.752 [-0.11, 0.11] 0.900 0.004 Very Strong Evidence for Null 1.000 7136 2748
is_A:pressure_self_cb -0.18 0.67 [-1.51, 1.16] 0.602 [-0.11, 0.11] 0.125 0.017 Very Strong Evidence for Null 1.001 994 1789
is_B:pressure_self_cb 0.11 0.57 [-1.02, 1.27] 0.581 [-0.11, 0.11] 0.168 0.017 Very Strong Evidence for Null 1.007 1017 1699
is_A:pressure_partner_cb 0.23 0.52 [-0.78, 1.28] 0.675 [-0.11, 0.11] 0.157 0.015 Very Strong Evidence for Null 1.003 1025 1644
is_B:pressure_partner_cb -0.12 0.73 [-1.59, 1.29] 0.570 [-0.11, 0.11] 0.127 0.016 Very Strong Evidence for Null 1.005 977 1542
is_A:day 0.41*** 0.07 [ 0.27, 0.55] 1.000 [-0.11, 0.11] 0.000 >100 Overwhelming Evidence 1.000 5190 2993
is_B:day 0.09 0.07 [-0.05, 0.23] 0.899 [-0.11, 0.11] 0.642 0.007 Very Strong Evidence for Null 1.001 5465 3359
sd(is_A) 0.71 0.09 [0.56, 0.93] NA NA NA NA NA 1.009 432 926
sd(is_B) 0.81 0.10 [0.63, 1.05] NA NA NA NA NA 1.005 610 888
sd(is_A:pushing_self_cw) 0.12 0.06 [0.01, 0.25] NA NA NA NA NA 1.002 814 588
sd(is_B:pushing_self_cw) 0.06 0.05 [0.00, 0.18] NA NA NA NA NA 1.003 1108 1004
sd(is_A:pushing_partner_cw) 0.13 0.05 [0.04, 0.24] NA NA NA NA NA 1.002 1571 1131
sd(is_B:pushing_partner_cw) 0.12 0.06 [0.02, 0.25] NA NA NA NA NA 1.002 1430 1166
sigma 0.87 0.01 [0.85, 0.89] NA NA NA NA NA 1.001 5267 3080

pushing

formula <- bf(
  aff ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:pushing_self_cw + is_B:pushing_self_cw + 
    is_A:pushing_partner_cw + is_B:pushing_partner_cw +
    
    is_A:pushing_self_cb + is_B:pushing_self_cb + 
    is_A:pushing_partner_cb + is_B:pushing_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    barriers_self_cw + barriers_self_cb + 
    is_A:day + is_B:day + 
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:pushing_self_cw + is_B:pushing_self_cw + 
    is_A:pushing_partner_cw + is_B:pushing_partner_cw |dd| coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
)


prior1 <- c(
  brms::set_prior("normal(0, 5)", class = "b")
  
  , brms::set_prior("normal(0, 20)", class = "b", coef = 'is_A')
  , brms::set_prior("normal(0, 20)", class = "b", coef = 'is_B')
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = gaussian()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

mood_gauss_pushing <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = gaussian(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("mood_gauss_pushing", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
if (check_models) {
  check_brms(mood_gauss_pushing, log_pp_check = FALSE)
  DHARMa.check_brms.all(mood_gauss_pushing, integer = FALSE)
}
## Model has no intercept. VIFs may not be sensible.
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                     Term  VIF   VIF 95% CI Increased SE Tolerance
##                     is_A 1.38 [1.32, 1.46]         1.18      0.72
##                     is_B 1.01 [1.00, 1.39]         1.01      0.99
##         barriers_self_cw 1.13 [1.09, 1.19]         1.06      0.89
##         barriers_self_cb 1.09 [1.06, 1.15]         1.05      0.91
##     is_A:pushing_self_cw 1.06 [1.03, 1.12]         1.03      0.95
##     is_B:pushing_self_cw 1.07 [1.04, 1.14]         1.04      0.93
##  is_A:pushing_partner_cw 1.12 [1.08, 1.18]         1.06      0.89
##  is_B:pushing_partner_cw 1.51 [1.44, 1.60]         1.23      0.66
##     is_A:pushing_self_cb 1.63 [1.55, 1.72]         1.28      0.61
##  Tolerance 95% CI
##      [0.68, 0.76]
##      [0.72, 1.00]
##      [0.84, 0.92]
##      [0.87, 0.95]
##      [0.89, 0.97]
##      [0.88, 0.96]
##      [0.85, 0.93]
##      [0.63, 0.70]
##      [0.58, 0.65]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        cauchy         41%
##        normal         41%
##   exponential          6%
## 
## Predicted Distribution of Response
## 
##   Distribution Probability
##        tweedie         41%
##  beta-binomial         38%
##    half-cauchy          6%
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## 
##  DHARMa outlier test based on exact binomial test with approximate
##  expectations
## 
## data:  model.check
## outliers at both margin(s) = 15, observations = 1776, p-value =
## 4.841e-06
## alternative hypothesis: true probability of success is not equal to 0.001998002
## 95 percent confidence interval:
##  0.004734615 0.013892099
## sample estimates:
## frequency of outliers (expected: 0.001998001998002 ) 
##                                          0.008445946
summary_mood <- summarize_brms(
  mood_gauss_pushing, 
  stats_to_report = stats_to_report,
  #model_rows_fixed = model_rows_fixed,
  #model_rows_random = model_rows_random,
  #model_rownames_fixed = model_rownames_fixed,
  #model_rownames_random = model_rownames_random,
  exponentiate = F) 
## Model has no intercept. VIFs may not be sensible.
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior
##   samples is recommended.
summary_mood %>%
  print_df(
    #rows_to_pack = rows_to_pack
    )
Est. SE 95% CI pd ROPE inside ROPE BF BF_Evidence Rhat Bulk_ESS Tail_ESS
is_A 3.65*** 0.12 [ 3.38, 3.89] 1.000 [-0.11, 0.11] 0.000 >100 Overwhelming Evidence 1.002 826 1559
is_B 3.84*** 0.15 [ 3.54, 4.13] 1.000 [-0.11, 0.11] 0.000 >100 Overwhelming Evidence 1.006 918 1668
barriers_self_cw -0.21*** 0.03 [-0.26, -0.16] 1.000 [-0.11, 0.11] 0.000 >100 Overwhelming Evidence 1.000 8770 3006
barriers_self_cb -0.35 0.22 [-0.80, 0.12] 0.936 [-0.11, 0.11] 0.122 0.050 Strong Evidence for Null 1.005 789 1517
is_A:pushing_self_cw 0.09 0.05 [-0.01, 0.18] 0.964 [-0.11, 0.11] 0.671 0.030 Strong Evidence for Null 1.002 4384 3110
is_B:pushing_self_cw -0.01 0.04 [-0.09, 0.06] 0.644 [-0.11, 0.11] 0.990 0.005 Very Strong Evidence for Null 1.000 5420 3240
is_A:pushing_partner_cw 0.09 0.05 [ 0.00, 0.18] 0.973 [-0.11, 0.11] 0.674 0.040 Strong Evidence for Null 1.000 3733 2884
is_B:pushing_partner_cw 0.12* 0.05 [ 0.01, 0.23] 0.984 [-0.11, 0.11] 0.387 0.077 Strong Evidence for Null 1.000 3298 2882
is_A:pushing_self_cb -0.17 0.49 [-1.19, 0.79] 0.637 [-0.11, 0.11] 0.165 0.014 Very Strong Evidence for Null 1.003 1189 1590
is_B:pushing_self_cb 0.82 0.68 [-0.53, 2.25] 0.891 [-0.11, 0.11] 0.065 0.043 Strong Evidence for Null 1.001 996 1543
is_A:pushing_partner_cb 0.78 0.56 [-0.35, 1.92] 0.916 [-0.11, 0.11] 0.054 0.040 Strong Evidence for Null 1.004 1020 1462
is_B:pushing_partner_cb -0.33 0.63 [-1.58, 1.01] 0.702 [-0.11, 0.11] 0.110 0.023 Very Strong Evidence for Null 1.002 1156 1592
is_A:day 0.35*** 0.10 [ 0.16, 0.54] 1.000 [-0.11, 0.11] 0.007 1.512 Weak Evidence 1.001 7719 2951
is_B:day 0.08 0.10 [-0.11, 0.28] 0.806 [-0.11, 0.11] 0.589 0.006 Very Strong Evidence for Null 1.004 6465 2807
sd(is_A) 0.66 0.09 [0.52, 0.87] NA NA NA NA NA 1.005 959 1597
sd(is_B) 0.80 0.10 [0.63, 1.05] NA NA NA NA NA 1.002 1041 1715
sd(is_A:pushing_self_cw) 0.10 0.07 [0.01, 0.24] NA NA NA NA NA 1.003 1278 1291
sd(is_B:pushing_self_cw) 0.07 0.06 [0.00, 0.20] NA NA NA NA NA 1.000 1483 1803
sd(is_A:pushing_partner_cw) 0.12 0.05 [0.02, 0.26] NA NA NA NA NA 1.001 1844 1640
sd(is_B:pushing_partner_cw) 0.13 0.07 [0.01, 0.30] NA NA NA NA NA 1.002 1520 1406
sigma 0.84 0.01 [0.81, 0.86] NA NA NA NA NA 1.001 5401 3024

Reactance

range(df_double$reactance, na.rm = T) 
## [1] 0 5
hist(df_double$reactance, breaks = 7) 

hist(log(df_double$reactance+0.1), breaks = 10)

introduce_binary_reactance <- function(data) {
  data$is_reactance <- factor(data$reactance > 0, levels = c(FALSE, TRUE), labels = c(0, 1))
  return(data)
}

df_double <- introduce_binary_reactance(df_double)
if (use_mi) {
  for (i in seq_along(implist)) {
    implist[[i]] <- introduce_binary_reactance(implist[[i]])
  }
}

Persuasion

formula <- bf(
  is_reactance ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:persuasion_self_cw + is_B:persuasion_self_cw + 
    is_A:persuasion_partner_cw + is_B:persuasion_partner_cw +
    
    is_A:persuasion_self_cb + is_B:persuasion_self_cb + 
    is_A:persuasion_partner_cb + is_B:persuasion_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    is_A:day + is_B:day + 
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:persuasion_self_cw + is_B:persuasion_self_cw + 
    is_A:persuasion_partner_cw + is_B:persuasion_partner_cw |dd| coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
  )



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")

  , brms::set_prior("normal(0, 10)", class = "b", coef = 'is_A')
  , brms::set_prior("normal(0, 10)", class = "b", coef = 'is_B')
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = bernoulli()
#  )



#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

is_reactance_persuasion <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::bernoulli(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("is_reactance_persuasion", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
if (check_models) {
  check_brms(is_reactance_persuasion)
  DHARMa.check_brms.all(is_reactance_persuasion, integer = FALSE)
}
## Model has no intercept. VIFs may not be sensible.
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                        Term  VIF   VIF 95% CI Increased SE Tolerance
##                        is_A 1.80 [1.69, 1.93]         1.34      0.55
##                        is_B 1.37 [1.30, 1.46]         1.17      0.73
##     is_A:persuasion_self_cw 1.17 [1.11, 1.24]         1.08      0.86
##     is_B:persuasion_self_cw 1.05 [1.02, 1.14]         1.02      0.95
##  is_A:persuasion_partner_cw 1.06 [1.02, 1.14]         1.03      0.94
##  is_B:persuasion_partner_cw 1.61 [1.51, 1.72]         1.27      0.62
##     is_A:persuasion_self_cb 2.00 [1.87, 2.15]         1.41      0.50
##  Tolerance 95% CI
##      [0.52, 0.59]
##      [0.68, 0.77]
##      [0.81, 0.90]
##      [0.88, 0.98]
##      [0.88, 0.98]
##      [0.58, 0.66]
##      [0.47, 0.54]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        normal         44%
##        cauchy         12%
##         gamma         12%
## 
## Predicted Distribution of Response
## 
##   Distribution Probability
##      bernoulli         81%
##  beta-binomial         12%
##       binomial          6%
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 9 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching
## for problematic observations.

## 
##  DHARMa outlier test based on exact binomial test with approximate
##  expectations
## 
## data:  model.check
## outliers at both margin(s) = 2, observations = 756, p-value =
## 0.6663
## alternative hypothesis: true probability of success is not equal to 0.001998002
## 95 percent confidence interval:
##  0.0003205434 0.0095234999
## sample estimates:
## frequency of outliers (expected: 0.001998001998002 ) 
##                                          0.002645503
summary_is_reactance <- summarize_brms(
  is_reactance_persuasion, 
  stats_to_report = stats_to_report,
  #model_rows_fixed = model_rows_fixed,
  #model_rows_random = model_rows_random,
  #model_rownames_fixed = model_rownames_fixed,
  #model_rownames_random = model_rownames_random,
  exponentiate = T) 
## Model has no intercept. VIFs may not be sensible.
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior
##   samples is recommended.
summary_is_reactance %>%
  print_df(
    #rows_to_pack = rows_to_pack
    )
OR SE 95% CI pd ROPE inside ROPE BF BF_Evidence Rhat Bulk_ESS Tail_ESS
is_A 0.18*** 0.09 [0.06, 0.42] 1.000 [0.83, 1.20] 0.000 14.783 Strong Evidence 1.000 2614 2734
is_B 0.29** 0.12 [0.12, 0.64] 0.999 [0.83, 1.20] 0.005 1.740 Weak Evidence 1.001 3369 3196
is_A:persuasion_self_cw 0.90 0.13 [0.64, 1.17] 0.769 [0.83, 1.20] 0.681 0.066 Strong Evidence for Null 1.001 2236 2875
is_B:persuasion_self_cw 0.78 0.11 [0.54, 1.02] 0.965 [0.83, 1.20] 0.308 0.309 Weak Evidence for Null 1.000 3027 2084
is_A:persuasion_partner_cw 0.88 0.17 [0.59, 1.28] 0.740 [0.83, 1.20] 0.573 0.084 Strong Evidence for Null 1.001 3281 2419
is_B:persuasion_partner_cw 1.26 0.24 [0.85, 1.93] 0.889 [0.83, 1.20] 0.370 0.154 Moderate Evidence for Null 1.001 2604 2514
is_A:persuasion_self_cb 5.18 4.43 [0.95, 36.77] 0.970 [0.83, 1.20] 0.029 0.520 Weak Evidence for Null 1.001 1970 2558
is_B:persuasion_self_cb 2.81 2.20 [0.58, 14.93] 0.905 [0.83, 1.20] 0.074 0.237 Moderate Evidence for Null 1.002 1937 2597
is_A:persuasion_partner_cb 1.54 1.10 [0.39, 8.13] 0.734 [0.83, 1.20] 0.175 0.106 Moderate Evidence for Null 1.002 1658 1952
is_B:persuasion_partner_cb 2.01 1.82 [0.33, 14.67] 0.784 [0.83, 1.20] 0.120 0.113 Moderate Evidence for Null 1.000 2156 2685
is_A:day 2.54 1.40 [0.88, 7.62] 0.959 [0.83, 1.20] 0.061 0.178 Moderate Evidence for Null 1.001 5826 3180
is_B:day 1.04 0.58 [0.34, 3.16] 0.528 [0.83, 1.20] 0.248 0.050 Strong Evidence for Null 1.003 5441 2933
sd(is_A) 1.37 0.42 [0.67, 2.43] NA NA NA NA NA 1.006 1215 1900
sd(is_B) 1.37 0.33 [0.84, 2.17] NA NA NA NA NA 1.000 2191 2654
sd(is_A:persuasion_self_cw) 0.30 0.18 [0.02, 0.71] NA NA NA NA NA 1.003 794 963
sd(is_B:persuasion_self_cw) 0.35 0.26 [0.02, 0.95] NA NA NA NA NA 1.002 797 1734
sd(is_A:persuasion_partner_cw) 0.48 0.27 [0.05, 1.16] NA NA NA NA NA 1.002 1413 1509
sd(is_B:persuasion_partner_cw) 0.64 0.28 [0.17, 1.37] NA NA NA NA NA 1.003 1425 1048

Pressure

formula <- bf(
  is_reactance ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:pressure_self_cw + is_B:pressure_self_cw + 
    is_A:pressure_partner_cw + is_B:pressure_partner_cw +
    
    is_A:pressure_self_cb + is_B:pressure_self_cb + 
    is_A:pressure_partner_cb + is_B:pressure_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    is_A:day + is_B:day + 
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:pressure_self_cw + is_B:pressure_self_cw + 
    is_A:pressure_partner_cw + is_B:pressure_partner_cw |dd| coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
  )



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")

  , brms::set_prior("normal(0, 10)", class = "b", coef = 'is_A')
  , brms::set_prior("normal(0, 10)", class = "b", coef = 'is_B')
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = bernoulli()
#  )



#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

is_reactance_pressure <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::bernoulli(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("is_reactance_pressure", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
if (check_models) {
  check_brms(is_reactance_pressure)
  DHARMa.check_brms.all(is_reactance_pressure, integer = FALSE)
}
## Model has no intercept. VIFs may not be sensible.
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                      Term  VIF   VIF 95% CI Increased SE Tolerance
##                      is_A 1.82 [1.70, 1.95]         1.35      0.55
##                      is_B 1.03 [1.01, 1.15]         1.02      0.97
##     is_A:pressure_self_cw 1.02 [1.00, 1.29]         1.01      0.98
##     is_B:pressure_self_cw 1.01 [1.00, 1.37]         1.01      0.99
##  is_A:pressure_partner_cw 1.01 [1.00, 2.53]         1.00      0.99
##  is_B:pressure_partner_cw 2.48 [2.30, 2.67]         1.57      0.40
##     is_A:pressure_self_cb 2.93 [2.71, 3.17]         1.71      0.34
##  Tolerance 95% CI
##      [0.51, 0.59]
##      [0.87, 0.99]
##      [0.78, 1.00]
##      [0.73, 1.00]
##      [0.39, 1.00]
##      [0.37, 0.43]
##      [0.32, 0.37]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        normal         34%
##          beta         16%
##        cauchy         12%
## 
## Predicted Distribution of Response
## 
##   Distribution Probability
##      bernoulli         81%
##  beta-binomial         12%
##       binomial          6%
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 26 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching
## for problematic observations.

## 
##  DHARMa outlier test based on exact binomial test with approximate
##  expectations
## 
## data:  model.check
## outliers at both margin(s) = 1, observations = 756, p-value = 1
## alternative hypothesis: true probability of success is not equal to 0.001998002
## 95 percent confidence interval:
##  0.0000334886 0.0073476538
## sample estimates:
## frequency of outliers (expected: 0.001998001998002 ) 
##                                          0.001322751
summary_is_reactance <- summarize_brms(
  is_reactance_pressure, 
  stats_to_report = stats_to_report,
  #model_rows_fixed = model_rows_fixed,
  #model_rows_random = model_rows_random,
  #model_rownames_fixed = model_rownames_fixed,
  #model_rownames_random = model_rownames_random,
  exponentiate = T) 
## Model has no intercept. VIFs may not be sensible.
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior
##   samples is recommended.
summary_is_reactance %>%
  print_df(
    #rows_to_pack = rows_to_pack
    )
OR SE 95% CI pd ROPE inside ROPE BF BF_Evidence Rhat Bulk_ESS Tail_ESS
is_A 0.15*** 0.06 [0.06, 0.32] 1.000 [0.83, 1.20] 0.000 87.472 Very Strong Evidence 1.001 2317 2099
is_B 0.21*** 0.07 [0.10, 0.39] 1.000 [0.83, 1.20] 0.000 95.261 Very Strong Evidence 1.000 3422 2844
is_A:pressure_self_cw 2.39 1.22 [0.79, 10.38] 0.951 [0.83, 1.20] 0.053 0.690 Weak Evidence for Null 1.001 1952 1596
is_B:pressure_self_cw 1.87 0.65 [0.82, 4.46] 0.943 [0.83, 1.20] 0.087 0.354 Weak Evidence for Null 1.000 2479 1818
is_A:pressure_partner_cw 0.87 0.66 [0.09, 3.64] 0.576 [0.83, 1.20] 0.202 0.127 Moderate Evidence for Null 1.002 2272 2073
is_B:pressure_partner_cw 1.55 0.76 [0.42, 5.01] 0.807 [0.83, 1.20] 0.176 0.129 Moderate Evidence for Null 1.001 3059 2828
is_A:pressure_self_cb 17.11 29.38 [0.80, 1025.83] 0.967 [0.83, 1.20] 0.015 0.601 Weak Evidence for Null 1.002 2273 1848
is_B:pressure_self_cb 22.64** 32.84 [1.90, 815.16] 0.996 [0.83, 1.20] 0.005 2.398 Weak Evidence 1.000 2655 2128
is_A:pressure_partner_cb 1.05 1.45 [0.05, 16.96] 0.514 [0.83, 1.20] 0.106 0.129 Moderate Evidence for Null 1.002 2319 1259
is_B:pressure_partner_cb 0.36 0.64 [0.00, 8.17] 0.725 [0.83, 1.20] 0.073 0.104 Moderate Evidence for Null 1.001 3171 2356
is_A:day 2.26 1.17 [0.81, 6.51] 0.943 [0.83, 1.20] 0.074 0.147 Moderate Evidence for Null 1.000 5731 2868
is_B:day 1.33 0.70 [0.47, 3.72] 0.705 [0.83, 1.20] 0.226 0.052 Strong Evidence for Null 1.002 4981 2411
sd(is_A) 1.46 0.39 [0.83, 2.40] NA NA NA NA NA 1.001 1456 2124
sd(is_B) 1.14 0.28 [0.66, 1.84] NA NA NA NA NA 1.000 1388 1435
sd(is_A:pressure_self_cw) 1.28 1.08 [0.05, 3.99] NA NA NA NA NA 1.003 872 1645
sd(is_B:pressure_self_cw) 1.02 0.57 [0.17, 2.51] NA NA NA NA NA 1.002 1315 1027
sd(is_A:pressure_partner_cw) 1.81 1.16 [0.15, 4.45] NA NA NA NA NA 1.001 1296 1321
sd(is_B:pressure_partner_cw) 0.71 0.64 [0.03, 2.84] NA NA NA NA NA 1.001 2176 1877

pushing

formula <- bf(
  is_reactance ~ 0 + 
    is_A + is_B + # Intercepts
    is_A:pushing_self_cw + is_B:pushing_self_cw + 
    is_A:pushing_partner_cw + is_B:pushing_partner_cw +
    
    is_A:pushing_self_cb + is_B:pushing_self_cb + 
    is_A:pushing_partner_cb + is_B:pushing_partner_cb +
    
    #is_A:plan_self + is_B:plan_self + 
    #is_A:plan_partner + is_B:plan_partner +
    barriers_self_cw + barriers_self_cb + 
    is_A:day + is_B:day + 
    
    # Random effects
    (0 + 
    is_A + is_B + # Intercepts
    is_A:pushing_self_cw + is_B:pushing_self_cw + 
    is_A:pushing_partner_cw + is_B:pushing_partner_cw |dd| coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
  )



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")

  , brms::set_prior("normal(0, 10)", class = "b", coef = 'is_A')
  , brms::set_prior("normal(0, 10)", class = "b", coef = 'is_B')
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = bernoulli()
#  )



#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

is_reactance_pushing <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::bernoulli(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("is_reactance_pushing", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
if (check_models) {
  check_brms(is_reactance_pushing)
  DHARMa.check_brms.all(is_reactance_pushing, integer = FALSE)
}
## Model has no intercept. VIFs may not be sensible.
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                     Term  VIF   VIF 95% CI Increased SE Tolerance
##                     is_A 1.89 [1.78, 2.01]         1.37      0.53
##                     is_B 1.07 [1.03, 1.13]         1.03      0.94
##         barriers_self_cw 1.11 [1.07, 1.17]         1.05      0.90
##         barriers_self_cb 1.24 [1.18, 1.31]         1.11      0.81
##     is_A:pushing_self_cw 1.12 [1.08, 1.18]         1.06      0.89
##     is_B:pushing_self_cw 1.02 [1.00, 1.21]         1.01      0.98
##  is_A:pushing_partner_cw 1.05 [1.02, 1.12]         1.03      0.95
##  is_B:pushing_partner_cw 1.47 [1.40, 1.56]         1.21      0.68
##     is_A:pushing_self_cb 1.37 [1.31, 1.45]         1.17      0.73
##  Tolerance 95% CI
##      [0.50, 0.56]
##      [0.88, 0.97]
##      [0.86, 0.94]
##      [0.77, 0.84]
##      [0.85, 0.93]
##      [0.82, 1.00]
##      [0.89, 0.98]
##      [0.64, 0.71]
##      [0.69, 0.77]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        normal         44%
##        cauchy         12%
##         gamma          9%
## 
## Predicted Distribution of Response
## 
##   Distribution Probability
##      bernoulli         75%
##  beta-binomial         16%
##       binomial          9%
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 8 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching
## for problematic observations.

## 
##  DHARMa outlier test based on exact binomial test with approximate
##  expectations
## 
## data:  model.check
## outliers at both margin(s) = 2, observations = 486, p-value =
## 0.2536
## alternative hypothesis: true probability of success is not equal to 0.001998002
## 95 percent confidence interval:
##  0.0004987621 0.0147859211
## sample estimates:
## frequency of outliers (expected: 0.001998001998002 ) 
##                                          0.004115226
summary_is_reactance <- summarize_brms(
  is_reactance_pushing, 
  stats_to_report = stats_to_report,
  #model_rows_fixed = model_rows_fixed,
  #model_rows_random = model_rows_random,
  #model_rownames_fixed = model_rownames_fixed,
  #model_rownames_random = model_rownames_random,
  exponentiate = T) 
## Model has no intercept. VIFs may not be sensible.
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior
##   samples is recommended.
summary_is_reactance %>%
  print_df(
    #rows_to_pack = rows_to_pack
    )
OR SE 95% CI pd ROPE inside ROPE BF BF_Evidence Rhat Bulk_ESS Tail_ESS
is_A 0.13*** 0.07 [0.04, 0.33] 1.000 [0.83, 1.20] 0.000 >100 Overwhelming Evidence 1.000 2447 2234
is_B 0.16*** 0.07 [0.06, 0.37] 1.000 [0.83, 1.20] 0.000 >100 Overwhelming Evidence 1.000 3196 2361
barriers_self_cw 0.78 0.14 [0.54, 1.10] 0.919 [0.83, 1.20] 0.342 0.139 Moderate Evidence for Null 1.000 5070 2919
barriers_self_cb 2.48 1.81 [0.62, 10.60] 0.900 [0.83, 1.20] 0.096 0.211 Moderate Evidence for Null 1.001 2695 2661
is_A:pushing_self_cw 1.35 0.25 [0.93, 2.02] 0.944 [0.83, 1.20] 0.249 0.232 Moderate Evidence for Null 1.000 4243 2725
is_B:pushing_self_cw 1.44 0.30 [0.97, 2.27] 0.964 [0.83, 1.20] 0.157 0.463 Weak Evidence for Null 1.001 3100 2454
is_A:pushing_partner_cw 0.77 0.24 [0.31, 1.33] 0.814 [0.83, 1.20] 0.342 0.119 Moderate Evidence for Null 1.003 2255 1751
is_B:pushing_partner_cw 1.16 0.24 [0.76, 1.74] 0.761 [0.83, 1.20] 0.500 0.080 Strong Evidence for Null 1.001 4260 3039
is_A:pushing_self_cb 29.91* 47.79 [1.80, 1723.06] 0.992 [0.83, 1.20] 0.008 1.715 Weak Evidence 1.002 1684 1694
is_B:pushing_self_cb 2.52 4.01 [0.11, 74.90] 0.725 [0.83, 1.20] 0.080 0.125 Moderate Evidence for Null 1.001 2147 2542
is_A:pushing_partner_cb 0.08 0.16 [0.00, 3.84] 0.898 [0.83, 1.20] 0.029 0.237 Moderate Evidence for Null 1.001 2349 2793
is_B:pushing_partner_cb 1.20 1.87 [0.04, 27.96] 0.545 [0.83, 1.20] 0.092 0.107 Moderate Evidence for Null 1.001 2505 2772
is_A:day 1.31 0.93 [0.32, 5.22] 0.658 [0.83, 1.20] 0.192 0.055 Strong Evidence for Null 1.002 4646 2904
is_B:day 1.38 0.92 [0.37, 5.25] 0.685 [0.83, 1.20] 0.199 0.065 Strong Evidence for Null 1.001 6392 2678
sd(is_A) 1.45 0.49 [0.61, 2.66] NA NA NA NA NA 1.000 1365 2100
sd(is_B) 1.30 0.40 [0.65, 2.26] NA NA NA NA NA 1.001 1836 2125
sd(is_A:pushing_self_cw) 0.25 0.22 [0.01, 0.86] NA NA NA NA NA 1.003 1342 1956
sd(is_B:pushing_self_cw) 0.57 0.31 [0.07, 1.39] NA NA NA NA NA 1.008 1046 1246
sd(is_A:pushing_partner_cw) 0.64 0.49 [0.04, 1.92] NA NA NA NA NA 1.005 1150 1424
sd(is_B:pushing_partner_cw) 0.29 0.25 [0.01, 1.06] NA NA NA NA NA 1.001 1769 2283 
report::report_system()

Analyses were conducted using the R Statistical language (version 4.4.2; R Core Team, 2024) on Windows 11 x64 (build 26100)

report::report_packages()
  • rmarkdown (version 2.29; Allaire J et al., 2024)
  • beepr (version 2.0; Bååth R, 2024)
  • R.methodsS3 (version 1.8.2; Bengtsson H, 2003)
  • R.oo (version 1.27.0; Bengtsson H, 2003)
  • R.utils (version 2.12.3; Bengtsson H, 2023)
  • brms (version 2.22.0; Bürkner P, 2017)
  • posterior (version 1.6.0; Bürkner P et al., 2024)
  • extrafont (version 0.19; Chang W, 2023)
  • digest (version 0.6.37; Eddelbuettel D, 2024)
  • Rcpp (version 1.0.13.1; Eddelbuettel D et al., 2024)
  • bayesplot (version 1.11.1; Gabry J, Mahr T, 2024)
  • lubridate (version 1.9.3; Grolemund G, Wickham H, 2011)
  • DHARMa (version 0.4.7; Hartig F, 2024)
  • rlang (version 1.1.4; Henry L, Wickham H, 2024)
  • tidybayes (version 3.0.7; Kay M, 2024)
  • wbCorr (version 0.1.22; Küng P, 2023)
  • see (version 0.9.0; Lüdecke D et al., 2021)
  • tibble (version 3.2.1; Müller K, Wickham H, 2023)
  • patchwork (version 1.3.0; Pedersen T, 2024)
  • R (version 4.4.2; R Core Team, 2024)
  • openxlsx (version 4.2.7.1; Schauberger P, Walker A, 2024)
  • plotly (version 4.10.4; Sievert C, 2020)
  • ggplot2 (version 3.5.1; Wickham H, 2016)
  • forcats (version 1.0.0; Wickham H, 2023)
  • stringr (version 1.5.1; Wickham H, 2023)
  • rvest (version 1.0.4; Wickham H, 2024)
  • tidyverse (version 2.0.0; Wickham H et al., 2019)
  • readxl (version 1.4.3; Wickham H, Bryan J, 2023)
  • dplyr (version 1.1.4; Wickham H et al., 2023)
  • purrr (version 1.0.2; Wickham H, Henry L, 2023)
  • readr (version 2.1.5; Wickham H et al., 2024)
  • xml2 (version 1.3.6; Wickham H et al., 2023)
  • tidyr (version 1.3.1; Wickham H et al., 2024)
  • ggridges (version 0.5.6; Wilke C, 2024)
  • knitr (version 1.49; Xie Y, 2024)
  • kableExtra (version 1.4.0; Zhu H, 2024)
report::cite_packages()